
> ##Replication File for "Authoritarian Audiences, Rhetoric, and Propaganda in International Crises: Evidence from China"
> ##Authors: Jessica Chen Weiss and Allan Dafoe
> 
> ##Replication Real History
> rm(list=ls())

> #Packages
> library(ggplot2)

> require(reshape2)

> require(coefplot)

> library(robustbase)

> setwd("~/Dropbox/Dafoe-Weiss-Empirics/18-04-12-ISQ-Raluca/Code/Replication File ISQ/")

> #------------------------------------------------------------------------------------------
> ##Producing dataset. 
> #------------------------------------------------------------------------------------------
> d <- read.csv("./MASTER CODE/150515RealHistorys.csv")

> #Drop first row
> orig.var.names <- d[1,]

> d <- d[-1,]

> #### Start and End Time
> d$start.time <- strptime(d$V8, "%Y-%m-%d %H:%M:%S")

> d$end.time <- strptime(d$V9, "%Y-%m-%d %H:%M:%S")

> #Drop those starting before Oct 15 2015 11pm
> d <- d[d$start.time>"2015-10-14 23:00:00",]

> ##Last start.time
> d$start.time[length(d[,1])]
[1] "2016-03-25 20:38:03 GMT"

> #Drop those who came in over the quota. All data is missing
> d <- d[d$quota!="overall",]

> ### Indicator for Partner
> #"Yes, those with an "rid" and "RISN" value are from the original panel (respondents should have both variables), and those with a "gid" and "sname" value are from the additional panel we brought on yesterday."
> d$partner <- rep(NA, length(d[,1]))

> d$partner[d$RISN=="" & d$gid!=""] <- "B"

> d$partner[d$RISN!="" & d$gid==""] <- "A"

> d$partner <- as.factor(d$partner)

> #Dropping those who didn't come from either partner
> which(is.na(d$partner)==TRUE)
 [1]   91   92   93  272  865  866  867  871 1095 1097 1099 1825 2815 3150 3294 3296 3358 4123 4131 4134 4198 5326 5509

> d<-d[!is.na(d$partner),]

> #### Attention Checks ####
> #Number of respondents who pass all attention checks.(only those who pass the first get to ctri)
> r <- d$term=="" & d$quota==""

> r2 <- (d$ctri== 2 | d$ctri==1) 

> table(r,r2)
       r2
r       FALSE TRUE
  FALSE  1151  460
  TRUE    622 3717

> #Number of respondents who answered asc
> 
> rc <- r==TRUE & r2==TRUE

> sum(rc)
[1] 3717

> d$att <- rep(0, length(d[,1]))

> d$att[rc] <- 1

> #Measure time taken
> d$time.taken <- difftime(d$end.time, d$start.time, units='mins')

> #### Creating Pre-Scenario Question Variables ####
> ## Pre Approval
> #Combining approval together
> 
> var <- d$as0

> var.r <- d$as0r

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> #summary(lm(out ~ order.r))
> d$asc0 <- out

> #d$asc0
> d$asc0.or <- order.r

> var <- d$ra0_2

> out <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> d$rac0 <- out

> var <- d$na1_1

> out <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> d$na1.v <- out

> var <- d$na2

> out <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> d$na2.v <- out

> #Recording variable so that about right is in the middle at 2
> d$na2.v[d$na2.v==3] <- -99

> d$na2.v[d$na2.v==2] <- 3

> d$na2.v[d$na2.v==-99] <- 2

> d$na2.v.o <- d$na2.v

> #Recording variable so that don't know and refuse to answer are also set to 2
> d$na2.v.dn <- rep(0, length(d[,1]))

> d$na2.v.dn[d$na2.v==8 | d$na2.v==9] <- 1

> d$na2.v.dn[d$na2.v==1 | d$na2.v==2 | d$na2.v==3] <- 0

> d$na2.v[d$na2.v==8] <- 2

> d$na2.v[d$na2.v==9] <- 2

> var <- d$na3

> out <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> d$na3.v <- out

> d$na3.v.dn <- rep(0, length(d[,1]))

> d$na3.v.dn[d$na3.v==8 | d$na3.v==9] <- 1

> d$na3.v.dn[d$na3.v!=8 & d$na3.v!=9] <- 0

> #moving 8 don't know; 9 refuse to answer to 3 moderate
> d$na3.v[d$na3.v==8] <- 3

> d$na3.v[d$na3.v==9] <- 3

> d$pre.questions <- rep(NA, length(d$na3.v))

> d$pre.questions[is.na(d$na3.v)] <- 0

> d$pre.questions[!is.na(d$na3.v)] <- 1

> #Last survey without a pre.question
> d$start.time[max(which(d$pre.questions==0))]
[1] "2016-03-14 15:56:57 GMT"

> #### Creating Main Variables ####
> 
> #### Creating Approval and Resolve Variables ####
> #d$as1: Disapprove first at the top
> #d$as1r: Approve first at the top
> 
> #Combining approval together
> var <- d$as1

> var.r <- d$as1r

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> d$asc <- out

> #d$asc
> #asc.or: asc is in order reversed
> d$asc.or <- order.r

> #Order effects for approval
> d$as1.n <- as.numeric(as.character(d$as1))

> d$as1.n2 <- as.numeric(as.character(d$as1r))

> #Measure Resolve
> #Combining resolve together
> var <- d$ra1_1

> var.r <- d$ra1r_1

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> d$rac <- out

> d$rac.or <- order.r

> #Note that anchoring doesn't seem to affect resolve answer
> d$ra1.n <- as.numeric(as.character(d$ra1_1))

> d$ra1.n2 <- as.numeric(as.character(d$ra1r_1))

> # Resolve 2
> d$ra3c <- as.numeric(as.character(d$ra3))

> #Perhaps this is reason to only use Resolve2 as a dichotomous variable (4 or not 4).
> d$ra3c.d <- d$ra3c

> d$ra3c.d[d$ra3c==4] <- 1

> d$ra3c.d[d$ra3c<4] <- 0

> #### Create treatment variables ####
> #his history
> # pro provocation
> # adiz adiz
> # adizp adiz adn defiance
> # eli-f elite cue framing
> 
> d$his.o <- d$his

> d$his <- as.numeric(as.character(d$his))

> r <- is.na(d$his)

> d$his[r] <- 0

> d$pro.o <- d$pro

> d$pro <- as.numeric(as.character(d$pro))

> r <- is.na(d$pro)

> d$pro[r] <- 0

> d$ADIZ.o <- d$ADIZ

> d$ADIZ <- as.numeric(as.character(d$ADIZ))

> r <- is.na(d$ADIZ)

> d$ADIZ[r] <- 0

> d$ADIZp.o <- d$ADIZp

> d$ADIZp <- as.numeric(as.character(d$ADIZp))

> r <- is.na(d$ADIZp)

> d$ADIZp[r] <- 0

> table(as.numeric(d$ADIZ), as.numeric(d$ADIZp))
   
       0    1
  0 3410 1278
  1 1262    0

> table(d$ADIZ, d$ADIZp)
   
       0    1
  0 3410 1278
  1 1262    0

> #This shows that d$ADIZ indicates ADIZ only, whereas d$ADIZp indicates only those who got both 
> #We want d$ADIZ to be 1 for d$ADIZp as well:
> d$ADIZ[d$ADIZp==1] <- 1

> table(d$ADIZ, d$ADIZp)
   
       0    1
  0 3410    0
  1 1262 1278

> d$eli.f.o <- d$eli.f

> d$eli.f <- as.numeric(as.character(d$eli.f))

> r <- is.na(d$eli.f)

> d$eli.f[r] <- 0

> d$eli.c.o <- d$eli.c

> d$eli.c <- as.numeric(as.character(d$eli.c))

> r <- is.na(d$eli.c)

> d$eli.c[r] <- 0

> median(d$na1.v, na.rm=T)
[1] 8

> d$na1.v2 <- d$na1.v

> d$na1.v2[is.na(d$na1.v)] <- 8 

> median(d$na2.v, na.rm=T)
[1] 2

> d$na2.v2 <- d$na2.v

> d$na2.v2[is.na(d$na2.v)] <- 3

> median(d$na3.v, na.rm=T)
[1] 3

> d$na3.v2 <- d$na3.v

> d$na3.v2[is.na(d$na3.v)] <- 3 

> median(d$asc0, na.rm=T)
[1] 4

> d$asc0.v2 <- d$asc0

> d$asc0.v2[is.na(d$asc0)] <- 4 

> median(d$rac0, na.rm=T)
[1] 51

> d$rac0.v2 <- d$rac0

> d$rac0.v2[is.na(d$rac0)] <- 50

> #Background
> #Gender, 1 is female, 0 is male
> d$gender <- as.numeric(as.character(d$d1))-1

> d$gender.o <- d$gender

> d$gender.m <- is.na(d$gender)

> d$gender[d$gender.m==TRUE] <- 1

> d$age <- 2015 - (as.numeric(as.character(d$d3))+1929)

> d$age.o <- d$age

> d$age.m <- is.na(d$age)

> d$age[d$age.m==TRUE] <- 30

> #"[d6] What is the highest level of education you have received?"
> # 01 No formal education 02 Elementary school 03 Middle school
> # 04 High school
> # 05 College 06 Masters 07 Doctoral
> var <- d$d6

> var.r <- d$d6r

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> d$educ <- out

> d$educ.o <- d$educ

> d$educ.m <- is.na(d$educ)

> d$educ[d$educ.m==TRUE] <- 3

> d$start.time.n <- as.numeric(d$start.time)-as.numeric(d$start.time[1])

> d$start.time.n2 <- (d$start.time.n)^2

> d$start.time.n3 <- (d$start.time.n)^3

> d$start.time.n4 <- (d$start.time.n)^4

> d$start.time.n5 <- (d$start.time.n)^5

> d$start.time.n6 <- (d$start.time.n)^6

> d$start.time.n7 <- (d$start.time.n)^7

> d$start.time.n8 <- (d$start.time.n)^8

> #Time dummy for after the first wave
> d$start.time.swd <- d$start.time > "2015-11-05 11:00:00 EST"

> d$st.o <- substr(d$start.time, 1, nchar(d$start.time)-4)

> d$st <- as.POSIXlt(d$st.o, tz="US/Mountain")

> ## cc Combining cc together
> var <- d$cc

> var.r <- d$ccr

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> d$ccc <- out

> d$cc.or <- order.r

> #Anchoring shifts answers by .26 out of 5
> summary(lm(d$ccc ~ d$cc.or))

Call:
lm(formula = d$ccc ~ d$cc.or)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.1460 -0.9008  0.0992  0.8540  2.0992 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  3.14600    0.02194 143.386  < 2e-16 ***
d$cc.or     -0.24520    0.03107  -7.893 3.72e-15 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.014 on 4262 degrees of freedom
  (1686 observations deleted due to missingness)
Multiple R-squared:  0.01441,	Adjusted R-squared:  0.01418 
F-statistic:  62.3 on 1 and 4262 DF,  p-value: 3.722e-15


> #Tests of Resolve
> #[tr] To what extent will other countries be more or less likely to challenge China in the future, if China doesn't defend its maritime interests more forcefully?
> var <- d$tr_1

> var.r <-d$trr_1

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> d$tr <- out

> d$tr.or <- order.r

> summary(lm(d$tr ~ d$tr.or))

Call:
lm(formula = d$tr ~ d$tr.or)

Residuals:
    Min      1Q  Median      3Q     Max 
-6.8844 -1.8677  0.1323  2.1156  3.1323 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  6.86766    0.05362  128.07   <2e-16 ***
d$tr.or      0.01672    0.07584    0.22    0.826    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.484 on 4289 degrees of freedom
  (1659 observations deleted due to missingness)
Multiple R-squared:  1.133e-05,	Adjusted R-squared:  -0.0002218 
F-statistic: 0.04861 on 1 and 4289 DF,  p-value: 0.8255


> hist(d$tr)

> hist(as.numeric(as.character(var)), col="blue")

> hist(as.numeric(as.character(var.r)), add=TRUE, alpha=.3)

> #[lf] To what extent have recent events in China’s surrounding waters made China lose face?
> d$lf.o <- d$lf

> d$lfr.o <- d$lfr

> var <- d$lf.o

> var.r <-d$lfr.o

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> d$lf <- out

> d$lf.or <- order.r

> summary(lm(d$lf ~ d$lf.or))

Call:
lm(formula = d$lf ~ d$lf.or)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.4757 -0.4757 -0.2261  0.5243  1.7739 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  2.22611    0.02019 110.245   <2e-16 ***
d$lf.or      0.24962    0.02857   8.738   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.9352 on 4285 degrees of freedom
  (1663 observations deleted due to missingness)
Multiple R-squared:  0.01751,	Adjusted R-squared:  0.01728 
F-statistic: 76.36 on 1 and 4285 DF,  p-value: < 2.2e-16


> hist(d$lf)

> hist(as.numeric(as.character(var)), col="blue")

> hist(as.numeric(as.character(var.r)), add=TRUE, alpha=.3)

> #Much more likely to select high levels of loss of face when high levels come first. 
> 
> 
> 
> # [nh] To what extent would you say that China’s national honor is at stake in this disagreement?
> # [0 national honor not at stake, 10 national honor very much at stake]
> 
> #[lf] To what extent have recent events in China’s surrounding waters made China lose face?
> d$nh.o <- d$nh_1

> d$nhr.o <- d$nhr_1

> var <- d$nh.o

> var.r <-d$nhr.o

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> d$nh <- out

> d$nh.or <- order.r

> summary(lm(d$nh ~ d$nh.or))

Call:
lm(formula = d$nh ~ d$nh.or)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.9490 -2.8868  0.1132  2.1132  5.1132 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  4.88681    0.06332  77.171   <2e-16 ***
d$nh.or      0.06221    0.08955   0.695    0.487    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.928 on 4274 degrees of freedom
  (1674 observations deleted due to missingness)
Multiple R-squared:  0.0001129,	Adjusted R-squared:  -0.0001211 
F-statistic: 0.4825 on 1 and 4274 DF,  p-value: 0.4873


> hist(d$nh)

> hist(as.numeric(as.character(var)), col="blue")

> hist(as.numeric(as.character(var.r)), add=TRUE, alpha=.3)

> #Lots select 0, and doesn't seem to be driven by anchor. 
> 
> 
> 
> #[thr] To what extent have recent events in China’s surrounding waters made China lose face?
> d$thr.o <- d$thr

> d$thrr.o <- d$thrr

> var <- d$thr.o

> var.r <-d$thrr.o

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> d$thr <- out

> d$thr.or <- order.r

> summary(lm(out ~ order.r))

Call:
lm(formula = out ~ order.r)

Residuals:
     Min       1Q   Median       3Q      Max 
-2.06133 -0.06133  0.08466  0.93867  1.08466 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  2.91534    0.01792  162.70  < 2e-16 ***
order.r      0.14599    0.02535    5.76 9.02e-09 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.8285 on 4272 degrees of freedom
  (1676 observations deleted due to missingness)
Multiple R-squared:  0.007705,	Adjusted R-squared:  0.007473 
F-statistic: 33.17 on 1 and 4272 DF,  p-value: 9.025e-09


> hist(out)

> hist(as.numeric(as.character(var)), col="blue")

> hist(as.numeric(as.character(var.r)), add=TRUE, alpha=.3)

> write.csv(d, "./MASTER CODE/real.csv")

> save(d, file="./MASTER CODE/real.Rdata")

> #------------------------------------------------------------------------------------------
> ## Constructing Publishable dataset. 
> #------------------------------------------------------------------------------------------
> 
> keepvars <- c("asc", "asc0", "pre.questions", "asc.or", "his", "pro", "ADIZ", "ADIZp", "eli.f",
+               "eli.c", "asc0.v2", "na1.v2", "na2.v2", "na3.v2", "na2.v.dn", "na3.v.dn", "start.time",
+               "start.time.n", "start.time.n2", "start.time.n3", "start.time.swd",
+               "na1.v", "na2.v", "na3.v", "att", "rac",  "rac.or", "ctri",  
+               "V1") ##variables that are not identifiable

> data.share <- d[, names(d) %in% keepvars] #keeping only shareable data

> write.csv(data.share, "./MASTER CODE/real_toshare.csv")

> d <- read.csv("./MASTER CODE/real_toshare.csv")

> #analysis data keeps only non-identifiable variables. 
> #generating random numbers of identifiable variables
> set.seed(2018)

> d$gender <- sample(c(0,1), size=nrow(d), replace=TRUE, prob=c(.37,.63))

> d$educ <- sample(seq(from = 1, to = 7, by = 1), size = nrow(d), replace = TRUE, prob=c(.0013, .0007, .417, .071, .452, .050, .007))

> d$age <- sample(c(11,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33 , 34,35,36,37,38,39,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,56,57,58,59,60,61,62,63,64,65,66,67,68,69, 70,72,73,76,80,82,85),
+                 size=nrow(d), replace=TRUE, prob=c(0.0002, 0.0005, 0.0052, 0.0087, 0.0118, 0.0124, 0.0153, 0.0148, 0.0156, 0.0319, 0.0207, 0.0266, 0.0284, 0.0240, 
+                                                    0.4316, 0.0237, 0.0242, 0.0286, 0.0205, 0.0328, 0.0182, 0.0183, 0.0148, 0.0158, 0.0163, 0.0136, 0.0126, 0.0126,
+                                                    0.0111, 0.0139, 0.0082, 0.0067, 0.0064, 0.0062, 0.0040, 0.0045, 0.0050, 0.0029, 0.0015, 0.0020, 0.0032, 0.0030,
+                                                    0.0024, 0.0025, 0.0020, 0.0018, 0.0012, 0.0027, 0.0024, 0.0015, 0.0015, 0.0007, 0.0007, 0.0003, 0.0007, 0.0002,
+                                                    0.0003, 0.0002, 0.0002, 0.0002, 0.0002))

> d$partner <- sample(c("A", "B"), size=nrow(d), replace=TRUE, prob=c(.295, .705)) 

> #Gen other variables 
> d$gender.o <- d$gender

> d$gender.m <- is.na(d$gender)

> d$gender[d$gender.m==TRUE] <- 1

> d$age.o <- d$age

> d$age.m <- is.na(d$age)

> d$age[d$age.m==TRUE] <- 30

> d$educ.o <- d$educ

> d$educ.m <- is.na(d$educ)

> d$educ[d$educ.m==TRUE] <- 3

> write.csv(d, "./MASTER CODE/realanalysis.csv")

> #-----------------------------------------------------------------------------------------
> #Replication Main Text and Appendix
> #-----------------------------------------------------------------------------------------
> {{format(Sys.Date(), format="%B %d %Y")}}
[1] "February 22 2019"

> #+ opts, include=FALSE,eval=TRUE
> knitr::opts_chunk$set(eval=TRUE, echo=TRUE, error=TRUE,message=FALSE, warning=FALSE,
+                       fig.width=15, fig.height=5,width=60,dev="pdf"
+ )

> rm(list = ls())

> setwd("~/Dropbox/Dafoe-Weiss-Empirics/18-04-12-ISQ-Raluca/Code/Replication File ISQ/")

> load("./MASTER CODE/real.Rdata")

> rb <- rep(1, length(d[,1]))

> #Figures location
> fig.loc <- "../../Figures/"

> #Packages
> require(ggplot2)

> require(reshape2)

> require(coefplot)

> require(robustbase)

> require(rms)

> require(sandwich)

> library(lubridate)

> library(stargazer)

> w <- 5

> h <- 4

> #finding value of Z that gives 10% one sided test
> xs <- seq(-2,0, 0.01)

> xs[which(pnorm(xs)-0.05==min(abs(pnorm(xs)-0.05)))]
[1] -1.64

> pnorm(-1.64)
[1] 0.05050258

> innerCI.n <- 1.64 #set so exclusion of CI implies two sided significance of 0.1

> outerCI.n <- 1.96 #set so exclusion of CI implies two sided significance of 0.05

> # Customized Stargazer Function, without asterisk inflation
> stargazerAD <- function(table, title="Default", filename="Figures/default.tex", dep.var.labels=NULL) {
+   a <- stargazer(table, title=title,  out=filename, 
+                  align=TRUE, star.cutoffs = c(0.1, 0.05, 0.01, 0.001), star.char = c("\\dagger", "*", "**", "***"), 
+                  notes="$^{\\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$", notes.append=FALSE, single.row=TRUE)
+   return(a)
+ }

> results.fun <- function(data, r, model, title.w, fig.loc, coefficients=NULL, innerCI.n=innerCI.n, outerCI.n=outerCI.n, w=10, h=6) {
+   d <- data
+   rob.results <- robcov(ols(model, dat=d[r,], x=TRUE))
+   m1 <-  lm(model,  dat=d[r,])
+   homo.results <- summary(m1)
+   n <- length(d[r,1])
+   title <- paste(title.w, ", n=", n, sep="")
+   pf <- coefplot(m1, title=title, 
+                  intercept=FALSE, lwdOuter=0.5, innerCI=innerCI.n, outerCI= outerCI.n,
+                  coefficients=coefficients) + theme_bw()
+   fig.results <- pf
+   list.results <- list("rob.results"=rob.results, "homo.results"=homo.results, "title"=title, "fig.results"=fig.results, "model"=m1, m1)
+   return(list.results)
+ }

> ## Number of observations ##
> n.obs <- sum(table(d$asc))

> cat(n.obs, file=paste(fig.loc,"Hist-n.tex", sep=""))

> n.obs.att <- sum(d$att) 

> cat(n.obs.att, file=paste(fig.loc,"Hist-n-att.tex", sep=""))

> #From PAP: 
> #For our primary specification we will estimate the effects of our manipulations using OLS with heteroscedasticity consis- tent standard errors, conditioning on the other experimental conditions to increase power. Hypothesis tests will be one-sided when appropriate. We will also report our estimates of the average causal effect when not conditioning on anything else, an estimand termed the average marginal component effect by Hainmueller, Hopkins, and Yamamoto.37
> #Our primary test of these hypotheses will involve OLS regression, controlling for the other experimental conditions, and in the real-history design for the pre-scenario answers to the same questions. 
> #robcov(ols(y ~ a.1 * a.2, toyData, x = TRUE)) ## heteroskedastic (ATE)
> 
> ### Results for Approval #### 
> #' ### Results for Approval #### 
> 
> 
> #r <-  rep(1, length(d[,1]))
> r <- !is.na(d$asc) & rb #& d$start.time > "2015-11-04 21:02:13 EST"

> type <- "Hist"

> outcome <- "Approval"

> cont <- ""

> title.w <- paste(type, ", Effect on ", outcome, cont, sep="")

> title.f <- paste(type, "-", outcome, cont, sep="")

> coefficients.list<- c("his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c")

> model <- asc ~ pre.questions + asc.or + his + pro + ADIZ + ADIZp + eli.f + eli.c

> results <- results.fun(data=d, r=r, model=model, title.w=title.w, fig.loc=fig.loc, coefficients=coefficients.list, innerCI.n=innerCI.n, outerCI.n=outerCI.n, w=w, h=h)

> results[1]
$rob.results
Linear Regression Model
 
 ols(formula = model, data = d[r, ], x = TRUE)
 
                Model Likelihood     Discrimination    
                   Ratio Test           Indexes        
 Obs    5445    LR chi2    143.37    R2       0.026    
 sigma1.0628    d.f.            8    R2 adj   0.025    
 d.f.   5436    Pr(> chi2) 0.0000    g        0.195    
 
 Residuals
 
     Min      1Q  Median      3Q     Max 
 -2.8432 -0.5743  0.3423  0.7458  1.8537 
 
 
               Coef    S.E.   t     Pr(>|t|)
 Intercept      3.3904 0.0464 73.09 <0.0001 
 pre.questions -0.1876 0.0436 -4.30 <0.0001 
 asc.or         0.3076 0.0288 10.67 <0.0001 
 his            0.0690 0.0334  2.06 0.0390  
 pro           -0.0565 0.0346 -1.63 0.1025  
 ADIZ           0.0514 0.0361  1.42 0.1544  
 ADIZp         -0.0352 0.0424 -0.83 0.4070  
 eli.f          0.0600 0.0551  1.09 0.2761  
 eli.c          0.0141 0.0561  0.25 0.8011  
 


> results[2] 
$homo.results

Call:
lm(formula = model, data = d[r, ])

Residuals:
    Min      1Q  Median      3Q     Max 
-2.8432 -0.5743  0.3423  0.7458  1.8537 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)    3.39043    0.04602  73.676  < 2e-16 ***
pre.questions -0.18761    0.04330  -4.333  1.5e-05 ***
asc.or         0.30760    0.02883  10.670  < 2e-16 ***
his            0.06903    0.03425   2.015   0.0439 *  
pro           -0.05650    0.03467  -1.629   0.1033    
ADIZ           0.05136    0.03620   1.419   0.1561    
ADIZp         -0.03518    0.04232  -0.831   0.4058    
eli.f          0.05999    0.05418   1.107   0.2682    
eli.c          0.01415    0.05548   0.255   0.7987    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.063 on 5436 degrees of freedom
Multiple R-squared:  0.02599,	Adjusted R-squared:  0.02455 
F-statistic: 18.13 on 8 and 5436 DF,  p-value: < 2.2e-16



> results$fig.results

> coefplot(lm(model, data=d[r,]), title=title.w, 
+          intercept=FALSE, lwdOuter=0.5, innerCI=innerCI.n, outerCI= outerCI.n, data=d) + theme_bw()

> m11 <- results[5]

> pdf(paste(fig.loc,title.f,".pdf", sep=""), width=w, height=h)

> results$fig.results

> dev.off()
RStudioGD 
        2 

> pro.results <- coef(summary(results$model))["pro",]

> pro.results <- rbind(pro.results, coef(summary(results$model))["ADIZp",])

> pro.results <- as.data.frame(pro.results)

> names(pro.results) <- c("est", "se", "t", "p")

> pro.results$var <- c("EP-3", "ADIZp")

> pro.results$dv <- "approval"

> pro.results$cov <- "no covariates"

> toprint <- paste("p_{+}=",round(pro.results[1,4]/2,3),  sep="")

> write(toprint, file=paste(fig.loc,type,outcome, "EP3", "p",".tex", sep=""))

> statement.results <- coef(summary(results$model))[c("ADIZ", "his", "eli.f", "eli.c"),]

> #Note: effect of ADIZ and ADIZp become positive after Nov 2015, consistent with any tough action being approved of.
> #Note, hypothesis tests should be read off of robcov function, to allow for heteroscedasticity. For convenience, the coefficient plot shows plots confidence intervals based on assuming homoscedasticity. 
> 
> 
> 
> 
> ### Now controlling for pre.questions, partner, time, and demographics ####
> #' ### Now controlling for pre.questions, partner, time, and demographics ###
> r <- !is.na(d$asc) & rb  #& d$start.time > "2015-11-04 21:02:13 EST"

> type <- "Hist"

> outcome <- "Approval"

> cont <- "Covariates"

> title.w <- paste(type, ", Effect on ", outcome, ", Controlling ", cont, sep="")

> title.f <- paste(type, "-", outcome, cont, sep="")

> coefficients.list<- c("his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c")

> model <- asc ~ his + pro + ADIZ + ADIZp + eli.f + eli.c +
+   partner + pre.questions + asc.or + asc0.v2  + na1.v2 + na2.v2 + na3.v2 + na2.v.dn + na3.v.dn + gender + educ + age + age.m + gender.m + educ.m + 
+   start.time.n + start.time.n2 + start.time.n3 + start.time.swd

> results <- results.fun(data=d, r=r, model=model, title.w=title.w, fig.loc=fig.loc, coefficients=coefficients.list, innerCI.n=innerCI.n, outerCI.n=outerCI.n, w=w, h=h)

> results[1]
$rob.results
Linear Regression Model
 
 ols(formula = model, data = d[r, ], x = TRUE)
 
                 Model Likelihood     Discrimination    
                    Ratio Test           Indexes        
 Obs    5445    LR chi2    1304.58    R2       0.213    
 sigma0.9567    d.f.            25    R2 adj   0.209    
 d.f.   5419    Pr(> chi2)  0.0000    g        0.544    
 
 Residuals
 
     Min      1Q  Median      3Q     Max 
 -3.1444 -0.5612  0.2081  0.6735  3.3489 
 
 
                Coef    S.E.         t    Pr(>|t|)
 Intercept       1.2585 1.003243e+20 0.00 1.0000  
 his             0.0717 4.191940e+18 0.00 1.0000  
 pro            -0.0542 1.201770e+18 0.00 1.0000  
 ADIZ            0.0665 7.376687e+17 0.00 1.0000  
 ADIZp          -0.0494 1.618164e+17 0.00 1.0000  
 eli.f           0.0878 1.254904e+18 0.00 1.0000  
 eli.c           0.0198 6.000028e+16 0.00 1.0000  
 partner=B      -0.0539 1.024188e+20 0.00 1.0000  
 pre.questions  -0.1675 4.515747e+19 0.00 1.0000  
 asc.or          0.3053 9.949937e+17 0.00 1.0000  
 asc0.v2         0.5363 3.588102e+18 0.00 1.0000  
 na1.v2         -0.0167 1.106783e+18 0.00 1.0000  
 na2.v2         -0.0820 2.311910e+18 0.00 1.0000  
 na3.v2          0.0133 4.848016e+17 0.00 1.0000  
 na2.v.dn        0.0417 1.768316e+18 0.00 1.0000  
 na3.v.dn       -0.0393 4.580653e+18 0.00 1.0000  
 gender         -0.0367 1.966822e+18 0.00 1.0000  
 educ            0.0414 1.886741e+18 0.00 1.0000  
 age             0.0027 9.474564e+16 0.00 1.0000  
 age.m           0.1368 2.545657e+33 0.00 1.0000  
 gender.m               6.803457e+19              
 educ.m          0.0736 2.396836e+14 0.00 1.0000  
 start.time.n    0.0000 1.756731e+07 0.00 1.0000  
 start.time.n2   0.0000 3.403000e-01 0.00 1.0000  
 start.time.n3   0.0000 6.769205e+20 0.00 1.0000  
 start.time.swd  0.0849 2.545657e+33 0.00 1.0000  
 


> results[2] 
$homo.results

Call:
lm(formula = model, data = d[r, ])

Residuals:
    Min      1Q  Median      3Q     Max 
-3.1444 -0.5612  0.2081  0.6735  3.3489 

Coefficients: (1 not defined because of singularities)
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)         1.258e+00  1.919e-01   6.559 5.91e-11 ***
his                 7.168e-02  3.088e-02   2.321 0.020327 *  
pro                -5.418e-02  3.124e-02  -1.735 0.082886 .  
ADIZ                6.647e-02  3.263e-02   2.037 0.041675 *  
ADIZp              -4.937e-02  3.812e-02  -1.295 0.195357    
eli.f               8.783e-02  4.882e-02   1.799 0.072028 .  
eli.c               1.976e-02  4.999e-02   0.395 0.692627    
partnerB           -5.389e-02  5.044e-02  -1.068 0.285363    
pre.questions      -1.675e-01  5.019e-02  -3.337 0.000852 ***
asc.or              3.053e-01  2.599e-02  11.748  < 2e-16 ***
asc0.v2             5.363e-01  1.558e-02  34.426  < 2e-16 ***
na1.v2             -1.668e-02  6.389e-03  -2.611 0.009061 ** 
na2.v2             -8.196e-02  2.166e-02  -3.784 0.000156 ***
na3.v2              1.328e-02  1.467e-02   0.905 0.365347    
na2.v.dn            4.165e-02  5.700e-02   0.731 0.464954    
na3.v.dn           -3.929e-02  7.899e-02  -0.497 0.618904    
gender             -3.668e-02  3.383e-02  -1.084 0.278295    
educ                4.138e-02  2.682e-02   1.543 0.122838    
age                 2.657e-03  1.640e-03   1.620 0.105316    
age.mTRUE           1.368e-01  4.804e-01   0.285 0.775862    
gender.mTRUE               NA         NA      NA       NA    
educ.mTRUE          7.360e-02  4.827e-01   0.152 0.878823    
start.time.n        1.274e-08  6.089e-08   0.209 0.834341    
start.time.n2      -5.227e-15  8.113e-15  -0.644 0.519405    
start.time.n3       3.412e-22  3.661e-22   0.932 0.351491    
start.time.swdTRUE  8.495e-02  1.214e-01   0.699 0.484286    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.9567 on 5420 degrees of freedom
Multiple R-squared:  0.2131,	Adjusted R-squared:  0.2096 
F-statistic: 61.14 on 24 and 5420 DF,  p-value: < 2.2e-16



> results$fig.results

> coefplot(lm(model, data=d[r,]), title=title.w, 
+          intercept=FALSE, lwdOuter=0.5, innerCI=innerCI.n, outerCI= outerCI.n, data=d) + theme_bw()

> m12 <- results[5]

> stargazerAD(c(m11, m12), title="Effect on Approval, History", filename=paste(fig.loc,title.f,"-table",".tex", sep=""))

% Table created by stargazer v.5.2.2 by Marek Hlavac, Harvard University. E-mail: hlavac at fas.harvard.edu
% Date and time: Fri, Feb 22, 2019 - 09:46:33
% Requires LaTeX packages: dcolumn 
\begin{table}[!htbp] \centering 
  \caption{Effect on Approval, History} 
  \label{} 
\begin{tabular}{@{\extracolsep{5pt}}lD{.}{.}{-3} D{.}{.}{-3} } 
\\[-1.8ex]\hline 
\hline \\[-1.8ex] 
 & \multicolumn{2}{c}{\textit{Dependent variable:}} \\ 
\cline{2-3} 
\\[-1.8ex] & \multicolumn{2}{c}{asc} \\ 
\\[-1.8ex] & \multicolumn{1}{c}{(1)} & \multicolumn{1}{c}{(2)}\\ 
\hline \\[-1.8ex] 
 pre.questions & -0.188^{***}$ $(0.043) & -0.167^{***}$ $(0.050) \\ 
  asc.or & 0.308^{***}$ $(0.029) & 0.305^{***}$ $(0.026) \\ 
  asc0.v2 &  & 0.536^{***}$ $(0.016) \\ 
  na1.v2 &  & -0.017^{**}$ $(0.006) \\ 
  na2.v2 &  & -0.082^{***}$ $(0.022) \\ 
  na3.v2 &  & 0.013$ $(0.015) \\ 
  na2.v.dn &  & 0.042$ $(0.057) \\ 
  na3.v.dn &  & -0.039$ $(0.079) \\ 
  gender &  & -0.037$ $(0.034) \\ 
  educ &  & 0.041$ $(0.027) \\ 
  age &  & 0.003$ $(0.002) \\ 
  age.m &  & 0.137$ $(0.480) \\ 
  gender.m &  &  \\ 
  educ.m &  & 0.074$ $(0.483) \\ 
  start.time.n &  & 0.000$ $(0.00000) \\ 
  start.time.n2 &  & -0.000$ $(0.000) \\ 
  start.time.n3 &  & 0.000$ $(0.000) \\ 
  start.time.swd &  & 0.085$ $(0.121) \\ 
  his & 0.069^{*}$ $(0.034) & 0.072^{*}$ $(0.031) \\ 
  pro & -0.056$ $(0.035) & -0.054^{\dagger}$ $(0.031) \\ 
  ADIZ & 0.051$ $(0.036) & 0.066^{*}$ $(0.033) \\ 
  ADIZp & -0.035$ $(0.042) & -0.049$ $(0.038) \\ 
  eli.f & 0.060$ $(0.054) & 0.088^{\dagger}$ $(0.049) \\ 
  eli.c & 0.014$ $(0.055) & 0.020$ $(0.050) \\ 
  partnerB &  & -0.054$ $(0.050) \\ 
  Constant & 3.390^{***}$ $(0.046) & 1.258^{***}$ $(0.192) \\ 
 \hline \\[-1.8ex] 
Observations & \multicolumn{1}{c}{5,445} & \multicolumn{1}{c}{5,445} \\ 
R$^{2}$ & \multicolumn{1}{c}{0.026} & \multicolumn{1}{c}{0.213} \\ 
Adjusted R$^{2}$ & \multicolumn{1}{c}{0.025} & \multicolumn{1}{c}{0.210} \\ 
Residual Std. Error & \multicolumn{1}{c}{1.063 (df = 5436)} & \multicolumn{1}{c}{0.957 (df = 5420)} \\ 
F Statistic & \multicolumn{1}{c}{18.129$^{***}$ (df = 8; 5436)} & \multicolumn{1}{c}{61.140$^{***}$ (df = 24; 5420)} \\ 
\hline 
\hline \\[-1.8ex] 
\textit{Note:}  & \multicolumn{2}{r}{$^{\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$} \\ 
\end{tabular} 
\end{table} 
 [1] ""                                                                                                                            
 [2] "% Table created by stargazer v.5.2.2 by Marek Hlavac, Harvard University. E-mail: hlavac at fas.harvard.edu"                 
 [3] "% Date and time: Fri, Feb 22, 2019 - 09:46:33"                                                                               
 [4] "% Requires LaTeX packages: dcolumn "                                                                                         
 [5] "\\begin{table}[!htbp] \\centering "                                                                                          
 [6] "  \\caption{Effect on Approval, History} "                                                                                   
 [7] "  \\label{} "                                                                                                                
 [8] "\\begin{tabular}{@{\\extracolsep{5pt}}lD{.}{.}{-3} D{.}{.}{-3} } "                                                           
 [9] "\\\\[-1.8ex]\\hline "                                                                                                        
[10] "\\hline \\\\[-1.8ex] "                                                                                                       
[11] " & \\multicolumn{2}{c}{\\textit{Dependent variable:}} \\\\ "                                                                 
[12] "\\cline{2-3} "                                                                                                               
[13] "\\\\[-1.8ex] & \\multicolumn{2}{c}{asc} \\\\ "                                                                               
[14] "\\\\[-1.8ex] & \\multicolumn{1}{c}{(1)} & \\multicolumn{1}{c}{(2)}\\\\ "                                                     
[15] "\\hline \\\\[-1.8ex] "                                                                                                       
[16] " pre.questions & -0.188^{***}$ $(0.043) & -0.167^{***}$ $(0.050) \\\\ "                                                      
[17] "  asc.or & 0.308^{***}$ $(0.029) & 0.305^{***}$ $(0.026) \\\\ "                                                              
[18] "  asc0.v2 &  & 0.536^{***}$ $(0.016) \\\\ "                                                                                  
[19] "  na1.v2 &  & -0.017^{**}$ $(0.006) \\\\ "                                                                                   
[20] "  na2.v2 &  & -0.082^{***}$ $(0.022) \\\\ "                                                                                  
[21] "  na3.v2 &  & 0.013$ $(0.015) \\\\ "                                                                                         
[22] "  na2.v.dn &  & 0.042$ $(0.057) \\\\ "                                                                                       
[23] "  na3.v.dn &  & -0.039$ $(0.079) \\\\ "                                                                                      
[24] "  gender &  & -0.037$ $(0.034) \\\\ "                                                                                        
[25] "  educ &  & 0.041$ $(0.027) \\\\ "                                                                                           
[26] "  age &  & 0.003$ $(0.002) \\\\ "                                                                                            
[27] "  age.m &  & 0.137$ $(0.480) \\\\ "                                                                                          
[28] "  gender.m &  &  \\\\ "                                                                                                      
[29] "  educ.m &  & 0.074$ $(0.483) \\\\ "                                                                                         
[30] "  start.time.n &  & 0.000$ $(0.00000) \\\\ "                                                                                 
[31] "  start.time.n2 &  & -0.000$ $(0.000) \\\\ "                                                                                 
[32] "  start.time.n3 &  & 0.000$ $(0.000) \\\\ "                                                                                  
[33] "  start.time.swd &  & 0.085$ $(0.121) \\\\ "                                                                                 
[34] "  his & 0.069^{*}$ $(0.034) & 0.072^{*}$ $(0.031) \\\\ "                                                                     
[35] "  pro & -0.056$ $(0.035) & -0.054^{\\dagger}$ $(0.031) \\\\ "                                                                
[36] "  ADIZ & 0.051$ $(0.036) & 0.066^{*}$ $(0.033) \\\\ "                                                                        
[37] "  ADIZp & -0.035$ $(0.042) & -0.049$ $(0.038) \\\\ "                                                                         
[38] "  eli.f & 0.060$ $(0.054) & 0.088^{\\dagger}$ $(0.049) \\\\ "                                                                
[39] "  eli.c & 0.014$ $(0.055) & 0.020$ $(0.050) \\\\ "                                                                           
[40] "  partnerB &  & -0.054$ $(0.050) \\\\ "                                                                                      
[41] "  Constant & 3.390^{***}$ $(0.046) & 1.258^{***}$ $(0.192) \\\\ "                                                            
[42] " \\hline \\\\[-1.8ex] "                                                                                                      
[43] "Observations & \\multicolumn{1}{c}{5,445} & \\multicolumn{1}{c}{5,445} \\\\ "                                                
[44] "R$^{2}$ & \\multicolumn{1}{c}{0.026} & \\multicolumn{1}{c}{0.213} \\\\ "                                                     
[45] "Adjusted R$^{2}$ & \\multicolumn{1}{c}{0.025} & \\multicolumn{1}{c}{0.210} \\\\ "                                            
[46] "Residual Std. Error & \\multicolumn{1}{c}{1.063 (df = 5436)} & \\multicolumn{1}{c}{0.957 (df = 5420)} \\\\ "                 
[47] "F Statistic & \\multicolumn{1}{c}{18.129$^{***}$ (df = 8; 5436)} & \\multicolumn{1}{c}{61.140$^{***}$ (df = 24; 5420)} \\\\ "
[48] "\\hline "                                                                                                                    
[49] "\\hline \\\\[-1.8ex] "                                                                                                       
[50] "\\textit{Note:}  & \\multicolumn{2}{r}{$^{\\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$} \\\\ "          
[51] "\\end{tabular} "                                                                                                             
[52] "\\end{table} "                                                                                                               

> smalld <- d[,c("asc", "asc0", "his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c", 
+                "partner", "asc.or", 
+                "na1.v", "na2.v", "na3.v", "na2.v.dn", "na3.v.dn", "gender.o", "educ.o", "age.o",
+                "age.m", "gender.m", "educ.m")]

> #, "start.time.n", "start.time.n2", "start.time.n3", "start.time.swd"
> stargazer(smalld, title="Summary Statistics", out=paste(fig.loc, "summary-statistics-history",".tex", sep=""))

% Table created by stargazer v.5.2.2 by Marek Hlavac, Harvard University. E-mail: hlavac at fas.harvard.edu
% Date and time: Fri, Feb 22, 2019 - 09:46:34
\begin{table}[!htbp] \centering 
  \caption{Summary Statistics} 
  \label{} 
\begin{tabular}{@{\extracolsep{5pt}}lccccccc} 
\\[-1.8ex]\hline 
\hline \\[-1.8ex] 
Statistic & \multicolumn{1}{c}{N} & \multicolumn{1}{c}{Mean} & \multicolumn{1}{c}{St. Dev.} & \multicolumn{1}{c}{Min} & \multicolumn{1}{c}{Pctl(25)} & \multicolumn{1}{c}{Pctl(75)} & \multicolumn{1}{c}{Max} \\ 
\hline \\[-1.8ex] 
asc & 5,445 & 3.405 & 1.076 & 1.000 & 3.000 & 4.000 & 5.000 \\ 
asc0 & 4,927 & 3.890 & 0.917 & 1.000 & 3.000 & 5.000 & 5.000 \\ 
his & 5,950 & 0.219 & 0.413 & 0 & 0 & 0 & 1 \\ 
pro & 5,950 & 0.211 & 0.408 & 0 & 0 & 0 & 1 \\ 
ADIZ & 5,950 & 0.427 & 0.495 & 0 & 0 & 1 & 1 \\ 
ADIZp & 5,950 & 0.215 & 0.411 & 0 & 0 & 0 & 1 \\ 
eli.f & 5,950 & 0.071 & 0.257 & 0 & 0 & 0 & 1 \\ 
eli.c & 5,950 & 0.067 & 0.251 & 0 & 0 & 0 & 1 \\ 
asc.or & 5,445 & 0.500 & 0.500 & 0.000 & 0.000 & 1.000 & 1.000 \\ 
na1.v & 4,791 & 7.698 & 2.275 & 0.000 & 6.000 & 10.000 & 10.000 \\ 
na2.v & 4,791 & 2.374 & 0.661 & 1.000 & 2.000 & 3.000 & 3.000 \\ 
na3.v & 4,791 & 3.143 & 0.958 & 1.000 & 3.000 & 4.000 & 5.000 \\ 
na2.v.dn & 5,950 & 0.065 & 0.247 & 0 & 0 & 0 & 1 \\ 
na3.v.dn & 5,950 & 0.032 & 0.176 & 0 & 0 & 0 & 1 \\ 
gender.o & 3,535 & 0.374 & 0.484 & 0.000 & 0.000 & 1.000 & 1.000 \\ 
educ.o & 3,531 & 4.941 & 0.606 & 1.000 & 5.000 & 5.000 & 7.000 \\ 
age.o & 3,535 & 34.331 & 10.166 & 11.000 & 27.000 & 40.000 & 85.000 \\ 
age.m & 5,950 & 0.406 & 0.491 & 0 & 0 & 1 & 1 \\ 
gender.m & 5,950 & 0.406 & 0.491 & 0 & 0 & 1 & 1 \\ 
educ.m & 5,950 & 0.407 & 0.491 & 0 & 0 & 1 & 1 \\ 
\hline \\[-1.8ex] 
\end{tabular} 
\end{table} 

> pdf(paste(fig.loc,title.f,".pdf", sep=""), width=w, height=h)

> results$fig.results

> dev.off()
RStudioGD 
        2 

> pro.results <- rbind(pro.results, 
+                      c(as.numeric(coef(summary(results$model))["pro",]), "EP-3", "approval", "covariates"))

> pro.results <- rbind(pro.results, 
+                      c(as.numeric(coef(summary(results$model))["ADIZp",]), "ADIZp", "approval", "covariates"))

> statement.results <- rbind(statement.results, coef(summary(results$model))[c("ADIZ", "his", "eli.f", "eli.c"),])

> row.names(statement.results) <- paste(row.names(statement.results), rep(", ", 4), c(rep("approval",4), rep("approval, cov", 4)), sep="")

> ##Analysis replicated among only those who passed attention check
> #r <-  rep(1, length(d[,1]))
> r <- !is.na(d$asc) & rb #& d$start.time > "2015-11-04 21:02:13 EST"

> type <- "Inattentive-Hist"

> outcome <- "Approval"

> cont <- ""

> title.w <- paste(type, ", Effect on ", outcome, cont, sep="")

> title.f <- paste(type, "-", outcome, cont, sep="")

> coefficients.list<- c("his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c")

> model <- asc ~ pre.questions + asc.or + his + pro + ADIZ + ADIZp + eli.f + eli.c

> results <- results.fun(data=d[d$att==1,], r=r, model=model, title.w=title.w, fig.loc=fig.loc, coefficients=coefficients.list, innerCI.n=innerCI.n, outerCI.n=outerCI.n, w=w, h=h)

> results[1]
$rob.results
Frequencies of Missing Values Due to Each Variable
          asc pre.questions        asc.or           his           pro          ADIZ         ADIZp         eli.f         eli.c 
         2016          2016          2016          2016          2016          2016          2016          2016          2016 

Linear Regression Model
 
 ols(formula = model, data = d[r, ], x = TRUE)
 
 
                Model Likelihood     Discrimination    
                   Ratio Test           Indexes        
 Obs    3429    LR chi2     74.27    R2       0.021    
 sigma1.0827    d.f.            8    R2 adj   0.019    
 d.f.   3420    Pr(> chi2) 0.0000    g        0.180    
 
 Residuals
 
     Min      1Q  Median      3Q     Max 
 -2.8005 -0.5903  0.3502  0.7738  1.8963 
 
 
               Coef    S.E.   t     Pr(>|t|)
 Intercept      3.4041 0.0586 58.06 <0.0001 
 pre.questions -0.2102 0.0555 -3.79 0.0002  
 asc.or         0.2627 0.0370  7.11 <0.0001 
 his            0.1013 0.0428  2.37 0.0179  
 pro           -0.0523 0.0444 -1.18 0.2387  
 ADIZ           0.0324 0.0466  0.70 0.4868  
 ADIZp         -0.0493 0.0544 -0.91 0.3640  
 eli.f          0.0410 0.0694  0.59 0.5543  
 eli.c         -0.0209 0.0727 -0.29 0.7740  
 


> results[2] 
$homo.results

Call:
lm(formula = model, data = d[r, ])

Residuals:
    Min      1Q  Median      3Q     Max 
-2.8005 -0.5903  0.3502  0.7738  1.8963 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)    3.40406    0.05854  58.153  < 2e-16 ***
pre.questions -0.21019    0.05489  -3.829 0.000131 ***
asc.or         0.26275    0.03701   7.099 1.52e-12 ***
his            0.10129    0.04415   2.294 0.021833 *  
pro           -0.05232    0.04425  -1.182 0.237114    
ADIZ           0.03237    0.04622   0.700 0.483766    
ADIZp         -0.04935    0.05386  -0.916 0.359643    
eli.f          0.04103    0.07010   0.585 0.558434    
eli.c         -0.02088    0.07157  -0.292 0.770494    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.083 on 3420 degrees of freedom
  (2016 observations deleted due to missingness)
Multiple R-squared:  0.02143,	Adjusted R-squared:  0.01914 
F-statistic:  9.36 on 8 and 3420 DF,  p-value: 7.331e-13



> results$fig.results

> m11_inattentive <- results[5]

> ##Analysis replicated among only those who passed attention check
> ### Now controlling for pre.questions, partner, time, and demographics ####
> #' ### Now controlling for pre.questions, partner, time, and demographics ###
> r <- !is.na(d$asc) & rb  #& d$start.time > "2015-11-04 21:02:13 EST"

> type <- "Inattentive-Hist"

> outcome <- "Approval"

> cont <- "Covariates"

> title.w <- paste(type, ", Effect on ", outcome, ", Controlling ", cont, sep="")

> title.f <- paste(type, "-", outcome, cont, sep="")

> coefficients.list<- c("his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c")

> model <- asc ~ his + pro + ADIZ + ADIZp + eli.f + eli.c +
+   partner + pre.questions + asc.or + asc0.v2  + na1.v2 + na2.v2 + na3.v2 + na2.v.dn + na3.v.dn + gender + educ + age + age.m + gender.m + educ.m + 
+   start.time.n + start.time.n2 + start.time.n3 + start.time.swd

> results <- results.fun(data=d[d$att==1,], r=r, model=model, title.w=title.w, fig.loc=fig.loc, coefficients=coefficients.list, innerCI.n=innerCI.n, outerCI.n=outerCI.n, w=w, h=h)

> results[1]
$rob.results
Frequencies of Missing Values Due to Each Variable
           asc            his            pro           ADIZ          ADIZp          eli.f          eli.c        partner 
          2016           2016           2016           2016           2016           2016           2016           2016 
 pre.questions         asc.or        asc0.v2         na1.v2         na2.v2         na3.v2       na2.v.dn       na3.v.dn 
          2016           2016           2016           2016           2016           2016           2016           2016 
        gender           educ            age          age.m       gender.m         educ.m   start.time.n  start.time.n2 
          2016           2016           2016           2016           2016           2016           2016           2016 
 start.time.n3 start.time.swd 
          2016           2016 

Linear Regression Model
 
 ols(formula = model, data = d[r, ], x = TRUE)
 
 
                Model Likelihood     Discrimination    
                   Ratio Test           Indexes        
 Obs    3429    LR chi2    964.95    R2       0.245    
 sigma0.9531    d.f.           25    R2 adj   0.240    
 d.f.   3403    Pr(> chi2) 0.0000    g        0.587    
 
 Residuals
 
     Min      1Q  Median      3Q     Max 
 -3.1804 -0.5546  0.2110  0.6587  3.4514 
 
 
                Coef    S.E.         t    Pr(>|t|)
 Intercept       1.1712 1.032762e+20 0.00 1.0000  
 his             0.0954 8.602009e+18 0.00 1.0000  
 pro            -0.0572 5.280774e+17 0.00 1.0000  
 ADIZ            0.0374 1.835553e+18 0.00 1.0000  
 ADIZp          -0.0750 1.602630e+18 0.00 1.0000  
 eli.f           0.1073 4.795293e+18 0.00 1.0000  
 eli.c          -0.0019 2.587048e+18 0.00 1.0000  
 partner=B      -0.0920 1.231743e+20 0.00 1.0000  
 pre.questions  -0.1798 5.942854e+19 0.00 1.0000  
 asc.or          0.2546 2.292346e+18 0.00 1.0000  
 asc0.v2         0.5917 5.913171e+18 0.00 1.0000  
 na1.v2         -0.0215 1.484927e+18 0.00 1.0000  
 na2.v2         -0.0889 4.762477e+18 0.00 1.0000  
 na3.v2          0.0033 7.733109e+17 0.00 1.0000  
 na2.v.dn        0.0298 1.531433e+18 0.00 1.0000  
 na3.v.dn       -0.0363 9.198821e+18 0.00 1.0000  
 gender         -0.0528 3.097674e+18 0.00 1.0000  
 educ            0.0393 1.583708e+18 0.00 1.0000  
 age             0.0029 1.438323e+17 0.00 1.0000  
 age.m           0.0876 9.961950e+33 0.00 1.0000  
 gender.m               9.158784e+19              
 educ.m          0.0930 2.742430e+14 0.00 1.0000  
 start.time.n    0.0000 1.783088e+07 0.00 1.0000  
 start.time.n2   0.0000 2.454000e-01 0.00 1.0000  
 start.time.n3   0.0000 8.405714e+20 0.00 1.0000  
 start.time.swd  0.0846 9.961950e+33 0.00 1.0000  
 


> results[2] 
$homo.results

Call:
lm(formula = model, data = d[r, ])

Residuals:
    Min      1Q  Median      3Q     Max 
-3.1804 -0.5546  0.2110  0.6587  3.4514 

Coefficients: (1 not defined because of singularities)
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)         1.171e+00  2.149e-01   5.451 5.37e-08 ***
his                 9.542e-02  3.898e-02   2.448  0.01441 *  
pro                -5.718e-02  3.900e-02  -1.466  0.14269    
ADIZ                3.742e-02  4.077e-02   0.918  0.35881    
ADIZp              -7.496e-02  4.747e-02  -1.579  0.11443    
eli.f               1.073e-01  6.192e-02   1.732  0.08331 .  
eli.c              -1.850e-03  6.312e-02  -0.029  0.97661    
partnerB           -9.203e-02  6.256e-02  -1.471  0.14136    
pre.questions      -1.798e-01  6.141e-02  -2.928  0.00343 ** 
asc.or              2.546e-01  3.266e-02   7.796 8.41e-15 ***
asc0.v2             5.917e-01  1.937e-02  30.552  < 2e-16 ***
na1.v2             -2.150e-02  8.268e-03  -2.601  0.00934 ** 
na2.v2             -8.888e-02  2.787e-02  -3.189  0.00144 ** 
na3.v2              3.268e-03  1.868e-02   0.175  0.86115    
na2.v.dn            2.980e-02  7.995e-02   0.373  0.70942    
na3.v.dn           -3.630e-02  1.197e-01  -0.303  0.76163    
gender             -5.281e-02  3.529e-02  -1.496  0.13462    
educ                3.931e-02  2.791e-02   1.408  0.15909    
age                 2.874e-03  1.723e-03   1.668  0.09539 .  
age.mTRUE           8.759e-02  4.839e-01   0.181  0.85639    
gender.mTRUE               NA         NA      NA       NA    
educ.mTRUE          9.300e-02  4.821e-01   0.193  0.84705    
start.time.n        4.851e-08  7.343e-08   0.661  0.50886    
start.time.n2      -1.140e-14  9.865e-15  -1.155  0.24797    
start.time.n3       6.689e-22  4.495e-22   1.488  0.13681    
start.time.swdTRUE  8.462e-02  1.529e-01   0.553  0.58004    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.9531 on 3404 degrees of freedom
  (2016 observations deleted due to missingness)
Multiple R-squared:  0.2453,	Adjusted R-squared:   0.24 
F-statistic:  46.1 on 24 and 3404 DF,  p-value: < 2.2e-16



> results$fig.results

> m12_inattentive <- results[5]

> stargazerAD(c(m11_inattentive, m12_inattentive), title="Effect on Approval, History (Inattentive respondents omitted)", filename=paste(fig.loc,title.f,"-table",".tex", sep=""))

% Table created by stargazer v.5.2.2 by Marek Hlavac, Harvard University. E-mail: hlavac at fas.harvard.edu
% Date and time: Fri, Feb 22, 2019 - 09:46:35
% Requires LaTeX packages: dcolumn 
\begin{table}[!htbp] \centering 
  \caption{Effect on Approval, History (Inattentive respondents omitted)} 
  \label{} 
\begin{tabular}{@{\extracolsep{5pt}}lD{.}{.}{-3} D{.}{.}{-3} } 
\\[-1.8ex]\hline 
\hline \\[-1.8ex] 
 & \multicolumn{2}{c}{\textit{Dependent variable:}} \\ 
\cline{2-3} 
\\[-1.8ex] & \multicolumn{2}{c}{asc} \\ 
\\[-1.8ex] & \multicolumn{1}{c}{(1)} & \multicolumn{1}{c}{(2)}\\ 
\hline \\[-1.8ex] 
 pre.questions & -0.210^{***}$ $(0.055) & -0.180^{**}$ $(0.061) \\ 
  asc.or & 0.263^{***}$ $(0.037) & 0.255^{***}$ $(0.033) \\ 
  asc0.v2 &  & 0.592^{***}$ $(0.019) \\ 
  na1.v2 &  & -0.022^{**}$ $(0.008) \\ 
  na2.v2 &  & -0.089^{**}$ $(0.028) \\ 
  na3.v2 &  & 0.003$ $(0.019) \\ 
  na2.v.dn &  & 0.030$ $(0.080) \\ 
  na3.v.dn &  & -0.036$ $(0.120) \\ 
  gender &  & -0.053$ $(0.035) \\ 
  educ &  & 0.039$ $(0.028) \\ 
  age &  & 0.003^{\dagger}$ $(0.002) \\ 
  age.m &  & 0.088$ $(0.484) \\ 
  gender.m &  &  \\ 
  educ.m &  & 0.093$ $(0.482) \\ 
  start.time.n &  & 0.00000$ $(0.00000) \\ 
  start.time.n2 &  & -0.000$ $(0.000) \\ 
  start.time.n3 &  & 0.000$ $(0.000) \\ 
  start.time.swd &  & 0.085$ $(0.153) \\ 
  his & 0.101^{*}$ $(0.044) & 0.095^{*}$ $(0.039) \\ 
  pro & -0.052$ $(0.044) & -0.057$ $(0.039) \\ 
  ADIZ & 0.032$ $(0.046) & 0.037$ $(0.041) \\ 
  ADIZp & -0.049$ $(0.054) & -0.075$ $(0.047) \\ 
  eli.f & 0.041$ $(0.070) & 0.107^{\dagger}$ $(0.062) \\ 
  eli.c & -0.021$ $(0.072) & -0.002$ $(0.063) \\ 
  partnerB &  & -0.092$ $(0.063) \\ 
  Constant & 3.404^{***}$ $(0.059) & 1.171^{***}$ $(0.215) \\ 
 \hline \\[-1.8ex] 
Observations & \multicolumn{1}{c}{3,429} & \multicolumn{1}{c}{3,429} \\ 
R$^{2}$ & \multicolumn{1}{c}{0.021} & \multicolumn{1}{c}{0.245} \\ 
Adjusted R$^{2}$ & \multicolumn{1}{c}{0.019} & \multicolumn{1}{c}{0.240} \\ 
Residual Std. Error & \multicolumn{1}{c}{1.083 (df = 3420)} & \multicolumn{1}{c}{0.953 (df = 3404)} \\ 
F Statistic & \multicolumn{1}{c}{9.360$^{***}$ (df = 8; 3420)} & \multicolumn{1}{c}{46.095$^{***}$ (df = 24; 3404)} \\ 
\hline 
\hline \\[-1.8ex] 
\textit{Note:}  & \multicolumn{2}{r}{$^{\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$} \\ 
\end{tabular} 
\end{table} 
 [1] ""                                                                                                                           
 [2] "% Table created by stargazer v.5.2.2 by Marek Hlavac, Harvard University. E-mail: hlavac at fas.harvard.edu"                
 [3] "% Date and time: Fri, Feb 22, 2019 - 09:46:35"                                                                              
 [4] "% Requires LaTeX packages: dcolumn "                                                                                        
 [5] "\\begin{table}[!htbp] \\centering "                                                                                         
 [6] "  \\caption{Effect on Approval, History (Inattentive respondents omitted)} "                                                
 [7] "  \\label{} "                                                                                                               
 [8] "\\begin{tabular}{@{\\extracolsep{5pt}}lD{.}{.}{-3} D{.}{.}{-3} } "                                                          
 [9] "\\\\[-1.8ex]\\hline "                                                                                                       
[10] "\\hline \\\\[-1.8ex] "                                                                                                      
[11] " & \\multicolumn{2}{c}{\\textit{Dependent variable:}} \\\\ "                                                                
[12] "\\cline{2-3} "                                                                                                              
[13] "\\\\[-1.8ex] & \\multicolumn{2}{c}{asc} \\\\ "                                                                              
[14] "\\\\[-1.8ex] & \\multicolumn{1}{c}{(1)} & \\multicolumn{1}{c}{(2)}\\\\ "                                                    
[15] "\\hline \\\\[-1.8ex] "                                                                                                      
[16] " pre.questions & -0.210^{***}$ $(0.055) & -0.180^{**}$ $(0.061) \\\\ "                                                      
[17] "  asc.or & 0.263^{***}$ $(0.037) & 0.255^{***}$ $(0.033) \\\\ "                                                             
[18] "  asc0.v2 &  & 0.592^{***}$ $(0.019) \\\\ "                                                                                 
[19] "  na1.v2 &  & -0.022^{**}$ $(0.008) \\\\ "                                                                                  
[20] "  na2.v2 &  & -0.089^{**}$ $(0.028) \\\\ "                                                                                  
[21] "  na3.v2 &  & 0.003$ $(0.019) \\\\ "                                                                                        
[22] "  na2.v.dn &  & 0.030$ $(0.080) \\\\ "                                                                                      
[23] "  na3.v.dn &  & -0.036$ $(0.120) \\\\ "                                                                                     
[24] "  gender &  & -0.053$ $(0.035) \\\\ "                                                                                       
[25] "  educ &  & 0.039$ $(0.028) \\\\ "                                                                                          
[26] "  age &  & 0.003^{\\dagger}$ $(0.002) \\\\ "                                                                                
[27] "  age.m &  & 0.088$ $(0.484) \\\\ "                                                                                         
[28] "  gender.m &  &  \\\\ "                                                                                                     
[29] "  educ.m &  & 0.093$ $(0.482) \\\\ "                                                                                        
[30] "  start.time.n &  & 0.00000$ $(0.00000) \\\\ "                                                                              
[31] "  start.time.n2 &  & -0.000$ $(0.000) \\\\ "                                                                                
[32] "  start.time.n3 &  & 0.000$ $(0.000) \\\\ "                                                                                 
[33] "  start.time.swd &  & 0.085$ $(0.153) \\\\ "                                                                                
[34] "  his & 0.101^{*}$ $(0.044) & 0.095^{*}$ $(0.039) \\\\ "                                                                    
[35] "  pro & -0.052$ $(0.044) & -0.057$ $(0.039) \\\\ "                                                                          
[36] "  ADIZ & 0.032$ $(0.046) & 0.037$ $(0.041) \\\\ "                                                                           
[37] "  ADIZp & -0.049$ $(0.054) & -0.075$ $(0.047) \\\\ "                                                                        
[38] "  eli.f & 0.041$ $(0.070) & 0.107^{\\dagger}$ $(0.062) \\\\ "                                                               
[39] "  eli.c & -0.021$ $(0.072) & -0.002$ $(0.063) \\\\ "                                                                        
[40] "  partnerB &  & -0.092$ $(0.063) \\\\ "                                                                                     
[41] "  Constant & 3.404^{***}$ $(0.059) & 1.171^{***}$ $(0.215) \\\\ "                                                           
[42] " \\hline \\\\[-1.8ex] "                                                                                                     
[43] "Observations & \\multicolumn{1}{c}{3,429} & \\multicolumn{1}{c}{3,429} \\\\ "                                               
[44] "R$^{2}$ & \\multicolumn{1}{c}{0.021} & \\multicolumn{1}{c}{0.245} \\\\ "                                                    
[45] "Adjusted R$^{2}$ & \\multicolumn{1}{c}{0.019} & \\multicolumn{1}{c}{0.240} \\\\ "                                           
[46] "Residual Std. Error & \\multicolumn{1}{c}{1.083 (df = 3420)} & \\multicolumn{1}{c}{0.953 (df = 3404)} \\\\ "                
[47] "F Statistic & \\multicolumn{1}{c}{9.360$^{***}$ (df = 8; 3420)} & \\multicolumn{1}{c}{46.095$^{***}$ (df = 24; 3404)} \\\\ "
[48] "\\hline "                                                                                                                   
[49] "\\hline \\\\[-1.8ex] "                                                                                                      
[50] "\\textit{Note:}  & \\multicolumn{2}{r}{$^{\\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$} \\\\ "         
[51] "\\end{tabular} "                                                                                                            
[52] "\\end{table} "                                                                                                              

> #analysis of distribution of inattentiveness 
> 
> 
> r <- !is.na(d$asc) & rb #& d$start.time > "2015-11-04 21:02:13 EST"

> type <- "Random-Inattentive-Hist"

> outcome <- "Approval"

> cont <- ""

> title.w <- paste(type, ", Effect on ", outcome, cont, sep="")

> title.f <- paste(type, "-", outcome, cont, sep="")

> coefficients.list<- c("his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c")

> model <- att ~ pre.questions + asc.or + his + pro + ADIZ + ADIZp + eli.f + eli.c

> results <- results.fun(data=d, r=r, model=model, title.w=title.w, fig.loc=fig.loc, coefficients=coefficients.list, innerCI.n=innerCI.n, outerCI.n=outerCI.n, w=w, h=h)

> results[1]
$rob.results
Linear Regression Model
 
 ols(formula = model, data = d[r, ], x = TRUE)
 
                Model Likelihood     Discrimination    
                   Ratio Test           Indexes        
 Obs    5445    LR chi2      6.50    R2       0.001    
 sigma0.4656    d.f.            8    R2 adj   0.000    
 d.f.   5436    Pr(> chi2) 0.5917    g        0.018    
 
 Residuals
 
     Min      1Q  Median      3Q     Max 
 -0.7182 -0.6663  0.3069  0.3212  0.3595 
 
 
               Coef    S.E.   t     Pr(>|t|)
 Intercept      0.6854 0.0201 34.07 <0.0001 
 pre.questions -0.0066 0.0189 -0.35 0.7245  
 asc.or        -0.0124 0.0126 -0.98 0.3248  
 his           -0.0085 0.0150 -0.56 0.5731  
 pro            0.0060 0.0151  0.40 0.6893  
 ADIZ           0.0268 0.0157  1.70 0.0885  
 ADIZp         -0.0010 0.0183 -0.06 0.9547  
 eli.f         -0.0174 0.0240 -0.73 0.4680  
 eli.c         -0.0131 0.0245 -0.53 0.5932  
 


> results[2] 
$homo.results

Call:
lm(formula = model, data = d[r, ])

Residuals:
    Min      1Q  Median      3Q     Max 
-0.7182 -0.6663  0.3069  0.3212  0.3595 

Coefficients:
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)    0.685399   0.020158  34.002   <2e-16 ***
pre.questions -0.006646   0.018965  -0.350    0.726    
asc.or        -0.012423   0.012628  -0.984    0.325    
his           -0.008452   0.015004  -0.563    0.573    
pro            0.006032   0.015188   0.397    0.691    
ADIZ           0.026805   0.015859   1.690    0.091 .  
ADIZp         -0.001039   0.018537  -0.056    0.955    
eli.f         -0.017401   0.023734  -0.733    0.463    
eli.c         -0.013067   0.024300  -0.538    0.591    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.4656 on 5436 degrees of freedom
Multiple R-squared:  0.001193,	Adjusted R-squared:  -0.0002773 
F-statistic: 0.8113 on 8 and 5436 DF,  p-value: 0.5925



> results$fig.results

> m31_inattentive <- results[5]

> ### Now controlling for pre.questions, partner, time, and demographics ####
> #' ### Now controlling for pre.questions, partner, time, and demographics ###
> r <- !is.na(d$asc) & rb  #& d$start.time > "2015-11-04 21:02:13 EST"

> type <- "Random-Inattentive-Hist"

> outcome <- "Approval"

> cont <- "Covariates"

> title.w <- paste(type, ", Effect on ", outcome, ", Controlling ", cont, sep="")

> title.f <- paste(type, "-", outcome, cont, sep="")

> coefficients.list<- c("his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c")

> model <- att ~ his + pro + ADIZ + ADIZp + eli.f + eli.c +
+   partner + pre.questions + asc.or + asc0.v2  + na1.v2 + na2.v2 + na3.v2 + na2.v.dn + na3.v.dn + gender + educ + age +
+   start.time.n + start.time.n2 + start.time.n3 + start.time.swd

> results <- results.fun(data=d, r=r, model=model, title.w=title.w, fig.loc=fig.loc, coefficients=coefficients.list, innerCI.n=innerCI.n, outerCI.n=outerCI.n, w=w, h=h)

> results[1]
$rob.results
Linear Regression Model
 
 ols(formula = model, data = d[r, ], x = TRUE)
 
                 Model Likelihood     Discrimination    
                    Ratio Test           Indexes        
 Obs    5445    LR chi2    6388.42    R2       0.691    
 sigma0.2594    d.f.            22    R2 adj   0.689    
 d.f.   5422    Pr(> chi2)  0.0000    g        0.429    
 
 Residuals
 
      Min       1Q   Median       3Q      Max 
 -1.38386 -0.19402 -0.02425  0.12536  1.46542 
 
 
                Coef    S.E.   t      Pr(>|t|)
 Intercept      -0.8831 0.0417 -21.19 <0.0001 
 his            -0.0090 0.0083  -1.08 0.2799  
 pro             0.0049 0.0084   0.58 0.5630  
 ADIZ            0.0108 0.0087   1.25 0.2120  
 ADIZp           0.0039 0.0102   0.38 0.7022  
 eli.f          -0.0096 0.0128  -0.76 0.4501  
 eli.c           0.0064 0.0142   0.45 0.6529  
 partner=B      -0.0157 0.0135  -1.16 0.2444  
 pre.questions   0.0223 0.0115   1.94 0.0527  
 asc.or         -0.0003 0.0070  -0.04 0.9690  
 asc0.v2         0.0000 0.0042   0.00 0.9980  
 na1.v2          0.0090 0.0020   4.47 <0.0001 
 na2.v2          0.0006 0.0059   0.10 0.9234  
 na3.v2          0.0143 0.0040   3.56 0.0004  
 na2.v.dn        0.0016 0.0188   0.09 0.9305  
 na3.v.dn       -0.0034 0.0281  -0.12 0.9032  
 gender         -0.1418 0.0088 -16.08 <0.0001 
 educ            0.3178 0.0057  55.69 <0.0001 
 age             0.0036 0.0005   8.01 <0.0001 
 start.time.n    0.0000 0.0000   2.26 0.0237  
 start.time.n2   0.0000 0.0000  -1.44 0.1489  
 start.time.n3   0.0000 0.0000   0.90 0.3664  
 start.time.swd -0.0440 0.0322  -1.36 0.1723  
 


> results[2] 
$homo.results

Call:
lm(formula = model, data = d[r, ])

Residuals:
     Min       1Q   Median       3Q      Max 
-1.38386 -0.19402 -0.02425  0.12536  1.46542 

Coefficients:
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)        -8.831e-01  3.954e-02 -22.335  < 2e-16 ***
his                -9.014e-03  8.373e-03  -1.077 0.281723    
pro                 4.856e-03  8.469e-03   0.573 0.566424    
ADIZ                1.082e-02  8.846e-03   1.223 0.221223    
ADIZp               3.915e-03  1.034e-02   0.379 0.704875    
eli.f              -9.634e-03  1.323e-02  -0.728 0.466654    
eli.c               6.402e-03  1.355e-02   0.472 0.636674    
partnerB           -1.570e-02  1.368e-02  -1.148 0.250929    
pre.questions       2.228e-02  1.361e-02   1.637 0.101669    
asc.or             -2.735e-04  7.046e-03  -0.039 0.969035    
asc0.v2            -1.049e-05  4.224e-03  -0.002 0.998020    
na1.v2              8.960e-03  1.730e-03   5.179 2.32e-07 ***
na2.v2              5.682e-04  5.872e-03   0.097 0.922925    
na3.v2              1.429e-02  3.976e-03   3.594 0.000328 ***
na2.v.dn            1.638e-03  1.544e-02   0.106 0.915481    
na3.v.dn           -3.421e-03  2.141e-02  -0.160 0.873046    
gender             -1.418e-01  8.778e-03 -16.153  < 2e-16 ***
educ                3.178e-01  4.064e-03  78.198  < 2e-16 ***
age                 3.626e-03  4.399e-04   8.243  < 2e-16 ***
start.time.n        3.360e-08  1.651e-08   2.035 0.041853 *  
start.time.n2      -3.021e-15  2.200e-15  -1.373 0.169746    
start.time.n3       9.014e-23  9.928e-23   0.908 0.363922    
start.time.swdTRUE -4.402e-02  3.292e-02  -1.337 0.181162    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.2594 on 5422 degrees of freedom
Multiple R-squared:  0.6906,	Adjusted R-squared:  0.6894 
F-statistic: 550.2 on 22 and 5422 DF,  p-value: < 2.2e-16



> results$fig.results

> m32_inattentive <- results[5]

> stargazerAD(c(m31_inattentive), title="Effect on Attentiveness, History", filename=paste(fig.loc,title.f,"-table",".tex", sep=""))

% Table created by stargazer v.5.2.2 by Marek Hlavac, Harvard University. E-mail: hlavac at fas.harvard.edu
% Date and time: Fri, Feb 22, 2019 - 09:46:37
% Requires LaTeX packages: dcolumn 
\begin{table}[!htbp] \centering 
  \caption{Effect on Attentiveness, History} 
  \label{} 
\begin{tabular}{@{\extracolsep{5pt}}lD{.}{.}{-3} } 
\\[-1.8ex]\hline 
\hline \\[-1.8ex] 
 & \multicolumn{1}{c}{\textit{Dependent variable:}} \\ 
\cline{2-2} 
\\[-1.8ex] & \multicolumn{1}{c}{att} \\ 
\hline \\[-1.8ex] 
 pre.questions & -0.007$ $(0.019) \\ 
  asc.or & -0.012$ $(0.013) \\ 
  his & -0.008$ $(0.015) \\ 
  pro & 0.006$ $(0.015) \\ 
  ADIZ & 0.027^{\dagger}$ $(0.016) \\ 
  ADIZp & -0.001$ $(0.019) \\ 
  eli.f & -0.017$ $(0.024) \\ 
  eli.c & -0.013$ $(0.024) \\ 
  Constant & 0.685^{***}$ $(0.020) \\ 
 \hline \\[-1.8ex] 
Observations & \multicolumn{1}{c}{5,445} \\ 
R$^{2}$ & \multicolumn{1}{c}{0.001} \\ 
Adjusted R$^{2}$ & \multicolumn{1}{c}{-0.0003} \\ 
Residual Std. Error & \multicolumn{1}{c}{0.466 (df = 5436)} \\ 
F Statistic & \multicolumn{1}{c}{0.811 (df = 8; 5436)} \\ 
\hline 
\hline \\[-1.8ex] 
\textit{Note:}  & \multicolumn{1}{r}{$^{\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$} \\ 
\end{tabular} 
\end{table} 
 [1] ""                                                                                                                  
 [2] "% Table created by stargazer v.5.2.2 by Marek Hlavac, Harvard University. E-mail: hlavac at fas.harvard.edu"       
 [3] "% Date and time: Fri, Feb 22, 2019 - 09:46:37"                                                                     
 [4] "% Requires LaTeX packages: dcolumn "                                                                               
 [5] "\\begin{table}[!htbp] \\centering "                                                                                
 [6] "  \\caption{Effect on Attentiveness, History} "                                                                    
 [7] "  \\label{} "                                                                                                      
 [8] "\\begin{tabular}{@{\\extracolsep{5pt}}lD{.}{.}{-3} } "                                                             
 [9] "\\\\[-1.8ex]\\hline "                                                                                              
[10] "\\hline \\\\[-1.8ex] "                                                                                             
[11] " & \\multicolumn{1}{c}{\\textit{Dependent variable:}} \\\\ "                                                       
[12] "\\cline{2-2} "                                                                                                     
[13] "\\\\[-1.8ex] & \\multicolumn{1}{c}{att} \\\\ "                                                                     
[14] "\\hline \\\\[-1.8ex] "                                                                                             
[15] " pre.questions & -0.007$ $(0.019) \\\\ "                                                                           
[16] "  asc.or & -0.012$ $(0.013) \\\\ "                                                                                 
[17] "  his & -0.008$ $(0.015) \\\\ "                                                                                    
[18] "  pro & 0.006$ $(0.015) \\\\ "                                                                                     
[19] "  ADIZ & 0.027^{\\dagger}$ $(0.016) \\\\ "                                                                         
[20] "  ADIZp & -0.001$ $(0.019) \\\\ "                                                                                  
[21] "  eli.f & -0.017$ $(0.024) \\\\ "                                                                                  
[22] "  eli.c & -0.013$ $(0.024) \\\\ "                                                                                  
[23] "  Constant & 0.685^{***}$ $(0.020) \\\\ "                                                                          
[24] " \\hline \\\\[-1.8ex] "                                                                                            
[25] "Observations & \\multicolumn{1}{c}{5,445} \\\\ "                                                                   
[26] "R$^{2}$ & \\multicolumn{1}{c}{0.001} \\\\ "                                                                        
[27] "Adjusted R$^{2}$ & \\multicolumn{1}{c}{-0.0003} \\\\ "                                                             
[28] "Residual Std. Error & \\multicolumn{1}{c}{0.466 (df = 5436)} \\\\ "                                                
[29] "F Statistic & \\multicolumn{1}{c}{0.811 (df = 8; 5436)} \\\\ "                                                     
[30] "\\hline "                                                                                                          
[31] "\\hline \\\\[-1.8ex] "                                                                                             
[32] "\\textit{Note:}  & \\multicolumn{1}{r}{$^{\\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$} \\\\ "
[33] "\\end{tabular} "                                                                                                   
[34] "\\end{table} "                                                                                                     

> ### Analyzing Resolve ####
> #' ### Analyzing Resolve ###
> r <- !is.na(d$rac) & rb #& d$start.time > "2015-11-04 21:02:13 EST"

> type <- "Hist"

> outcome <- "Resolve"

> cont <- ""

> title.w <- paste(type, ", Effect on ", outcome, cont, sep="")

> title.f <- paste(type, "-", outcome, cont, sep="")

> coefficients.list<- c("his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c")

> model <- rac ~ pre.questions + rac.or + his + pro + ADIZ + ADIZp + eli.f + eli.c

> results <- results.fun(data=d, r=r, model=model, title.w=title.w, fig.loc=fig.loc, coefficients=coefficients.list, innerCI.n=innerCI.n, outerCI.n=outerCI.n, w=w, h=h)

> results[1]
$rob.results
Linear Regression Model
 
 ols(formula = model, data = d[r, ], x = TRUE)
 
                 Model Likelihood     Discrimination    
                    Ratio Test           Indexes        
 Obs     5287    LR chi2     15.63    R2       0.003    
 sigma24.4907    d.f.            8    R2 adj   0.001    
 d.f.    5278    Pr(> chi2) 0.0479    g        1.486    
 
 Residuals
 
      Min       1Q   Median       3Q      Max 
 -55.5232 -15.5506  -0.8213  17.7325  50.7229 
 
 
               Coef    S.E.   t     Pr(>|t|)
 Intercept     50.7292 1.0616 47.79 <0.0001 
 pre.questions  2.0412 0.9913  2.06 0.0395  
 rac.or         1.2622 0.6733  1.87 0.0609  
 his           -0.9492 0.8117 -1.17 0.2423  
 pro            1.4906 0.8102  1.84 0.0659  
 ADIZ          -0.5030 0.8466 -0.59 0.5525  
 ADIZp          1.5005 0.9955  1.51 0.1318  
 eli.f         -0.7610 1.2536 -0.61 0.5439  
 eli.c         -0.8697 1.3724 -0.63 0.5263  
 


> results[2]
$homo.results

Call:
lm(formula = model, data = d[r, ])

Residuals:
    Min      1Q  Median      3Q     Max 
-55.523 -15.551  -0.821  17.733  50.723 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)    50.7292     1.0806  46.947   <2e-16 ***
pre.questions   2.0412     1.0162   2.009   0.0446 *  
rac.or          1.2622     0.6738   1.873   0.0611 .  
his            -0.9492     0.8001  -1.186   0.2355    
pro             1.4906     0.8095   1.841   0.0656 .  
ADIZ           -0.5030     0.8470  -0.594   0.5527    
ADIZp           1.5005     0.9872   1.520   0.1286    
eli.f          -0.7610     1.2612  -0.603   0.5463    
eli.c          -0.8697     1.2981  -0.670   0.5029    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 24.49 on 5278 degrees of freedom
Multiple R-squared:  0.002953,	Adjusted R-squared:  0.001441 
F-statistic: 1.954 on 8 and 5278 DF,  p-value: 0.04823



> results$fig.results

> coefplot(lm(model, data=d[r,]), title=title.w, 
+          intercept=FALSE, lwdOuter=0.5, innerCI=innerCI.n, outerCI= outerCI.n, data=d) + theme_bw()

> pdf(paste(fig.loc,title.f,".pdf", sep=""), width=w, height=h)

> results$fig.results

> dev.off()
RStudioGD 
        2 

> pro.results <- rbind(pro.results, 
+                      c(as.numeric(coef(summary(results$model))["pro",]), "EP-3", "resolve", "no covariates"))

> pro.results <- rbind(pro.results, 
+                      c(as.numeric(coef(summary(results$model))["ADIZp",]), "ADIZp", "resolve", "no covariates"))

> toprint <- paste("p_{+}=",round(as.numeric(pro.results[5,4])/2,3),  sep="")

> write(toprint, file=paste(fig.loc,type,outcome, "EP3", "p",".tex", sep=""))

> toprint <- paste("p_{+}=",round(as.numeric(pro.results[6,4])/2,3),  sep="")

> write(toprint, file=paste(fig.loc,type,outcome, "ADIZp", "p",".tex", sep=""))

> statement.results <- rbind(statement.results, coef(summary(results$model))[c("ADIZ", "his", "eli.f", "eli.c"),])

> r <- 9:12

> row.names(statement.results)[r] <- paste(row.names(statement.results)[r], rep(", ", 4), c(rep("resolve",4)), sep="")

> r <- !is.na(d$rac) & rb #& d$start.time > "2015-11-04 21:02:13 EST"

> type <- "Hist"

> outcome <- "Resolve"

> cont <- "Covariates"

> title.w <- paste(type, ", Effect on ", outcome, ", Controlling ", cont, sep="")

> title.f <- paste(type, "-", outcome, cont, sep="")

> coefficients.list<- c("his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c")

> model <- rac ~ pre.questions + his + pro + ADIZ + ADIZp + eli.f + eli.c +
+   partner + pre.questions + rac0.v2 + rac.or + na1.v2 + na2.v2 + na3.v2 + na2.v.dn + na3.v.dn +
+   gender + educ + age + age.m + gender.m + educ.m + 
+   start.time.n + start.time.n2 + start.time.n3 + start.time.swd

> results <- results.fun(data=d, r=r, model=model, title.w=title.w, fig.loc=fig.loc, coefficients=coefficients.list, innerCI.n=innerCI.n, outerCI.n=outerCI.n, w=w, h=h)

> results[1]
$rob.results
Linear Regression Model
 
 ols(formula = model, data = d[r, ], x = TRUE)
 
                  Model Likelihood     Discrimination    
                     Ratio Test           Indexes        
 Obs     5287    LR chi2    4456.36    R2       0.570    
 sigma16.1166    d.f.            25    R2 adj   0.567    
 d.f.    5261    Pr(> chi2)  0.0000    g       20.769    
 
 Residuals
 
     Min      1Q  Median      3Q     Max 
 -93.232  -7.306  -0.564   7.157  90.907 
 
 
                Coef    S.E.         t    Pr(>|t|)
 Intercept       3.8946 2.024104e+21 0.00 1.0000  
 pre.questions  -0.0032 7.581626e+20 0.00 1.0000  
 his            -0.2321 7.438555e+19 0.00 1.0000  
 pro             0.6813 1.877211e+19 0.00 1.0000  
 ADIZ           -1.0097 1.700121e+19 0.00 1.0000  
 ADIZp           0.6836 2.128897e+18 0.00 1.0000  
 eli.f           0.0411 1.521669e+19 0.00 1.0000  
 eli.c           0.6698 2.226933e+19 0.00 1.0000  
 partner=B       1.1264 1.672818e+21 0.00 1.0000  
 rac0.v2         0.7880 1.938319e+17 0.00 1.0000  
 rac.or          0.8427 2.518997e+17 0.00 1.0000  
 na1.v2          0.7214 1.176622e+19 0.00 1.0000  
 na2.v2         -0.2812 4.395259e+19 0.00 1.0000  
 na3.v2         -0.1528 6.298064e+18 0.00 1.0000  
 na2.v.dn       -0.4376 5.616407e+19 0.00 1.0000  
 na3.v.dn       -2.4840 9.989136e+19 0.00 1.0000  
 gender          0.9837 7.773077e+19 0.00 1.0000  
 educ            0.3434 8.849096e+19 0.00 1.0000  
 age             0.0269 4.865554e+16 0.00 1.0000  
 age.m           1.7289 1.283271e+35 0.00 1.0000  
 gender.m               1.117333e+21              
 educ.m         -0.9989 3.902068e+15 0.00 1.0000  
 start.time.n    0.0000 2.846023e+08 0.00 1.0000  
 start.time.n2   0.0000 5.453600e+00 0.00 1.0000  
 start.time.n3   0.0000 1.106882e+22 0.00 1.0000  
 start.time.swd -2.1125 1.283271e+35 0.00 1.0000  
 


> results[2]
$homo.results

Call:
lm(formula = model, data = d[r, ])

Residuals:
    Min      1Q  Median      3Q     Max 
-93.232  -7.306  -0.564   7.157  90.907 

Coefficients: (1 not defined because of singularities)
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)         3.895e+00  3.131e+00   1.244   0.2136    
pre.questions      -3.156e-03  8.582e-01  -0.004   0.9971    
his                -2.321e-01  5.276e-01  -0.440   0.6600    
pro                 6.813e-01  5.332e-01   1.278   0.2014    
ADIZ               -1.010e+00  5.580e-01  -1.809   0.0704 .  
ADIZp               6.836e-01  6.502e-01   1.051   0.2931    
eli.f               4.115e-02  8.306e-01   0.050   0.9605    
eli.c               6.698e-01  8.554e-01   0.783   0.4336    
partnerB            1.126e+00  8.597e-01   1.310   0.1902    
rac0.v2             7.880e-01  1.004e-02  78.461  < 2e-16 ***
rac.or              8.427e-01  4.438e-01   1.899   0.0576 .  
na1.v2              7.214e-01  1.097e-01   6.577 5.25e-11 ***
na2.v2             -2.812e-01  3.696e-01  -0.761   0.4468    
na3.v2             -1.528e-01  2.503e-01  -0.611   0.5415    
na2.v.dn           -4.376e-01  9.810e-01  -0.446   0.6556    
na3.v.dn           -2.484e+00  1.358e+00  -1.829   0.0675 .  
gender              9.837e-01  5.697e-01   1.727   0.0843 .  
educ                3.434e-01  4.520e-01   0.760   0.4475    
age                 2.689e-02  2.762e-02   0.974   0.3302    
age.mTRUE           1.729e+00  8.095e+00   0.214   0.8309    
gender.mTRUE               NA         NA      NA       NA    
educ.mTRUE         -9.989e-01  8.132e+00  -0.123   0.9022    
start.time.n        3.748e-07  1.037e-06   0.361   0.7178    
start.time.n2      -1.121e-14  1.383e-13  -0.081   0.9354    
start.time.n3      -1.171e-22  6.244e-21  -0.019   0.9850    
start.time.swdTRUE -2.112e+00  2.071e+00  -1.020   0.3078    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 16.12 on 5262 degrees of freedom
Multiple R-squared:  0.5695,	Adjusted R-squared:  0.5676 
F-statistic: 290.1 on 24 and 5262 DF,  p-value: < 2.2e-16



> results$fig.results

> coefplot(lm(model, data=d[r,]), title=title.w, 
+          intercept=FALSE, lwdOuter=0.5, innerCI=innerCI.n, outerCI= outerCI.n, data=d) + theme_bw()

> pdf(paste(fig.loc,title.f,".pdf", sep=""), width=w, height=h)

> results$fig.results

> dev.off()
RStudioGD 
        2 

> pro.results <- rbind(pro.results, 
+                      c(as.numeric(coef(summary(results$model))["pro",]), "EP-3", "resolve", "covariates"))

> pro.results <- rbind(pro.results, 
+                      c(as.numeric(coef(summary(results$model))["ADIZp",]), "ADIZp", "resolve", "covariates"))

> toprint <- paste("p_{+}=",round(as.numeric(pro.results[7,4])/2,3),  sep="")

> write(toprint, file=paste(fig.loc,type,outcome, "EP3", "p", "Cov", ".tex", sep=""))

> toprint <- paste("p_{+}=",round(as.numeric(pro.results[8,4])/2,3),  sep="")

> write(toprint, file=paste(fig.loc,type,outcome, "ADIZp", "p","Cov",".tex", sep=""))

> statement.results <- rbind(statement.results, coef(summary(results$model))[c("ADIZ", "his", "eli.f", "eli.c"),])

> r <- 13:16

> row.names(statement.results)[r] <- paste(row.names(statement.results)[r], rep(", ", 4), c(rep("resolve, cov",4)), sep="")

> ## Distribution of Approval ##
> 
> hist(d$asc0, main="Historical Approval (Pre)", breaks=seq(0.5, 5.5, 1), freq=FALSE, labels=TRUE, col="blue", xaxt='n', xlab="", cex=0.8)

> axis(side=1, at=c(1, 3, 5), labels=c("Strongly \n disapprove", "Neither approve \n nor disapprove", "Strongly \n approve"))

> pdf(paste(fig.loc,"Hist-Approval(Pre)-Histogram",".pdf", sep=""), width=w, height=h*1.6)

> par(mar=c(3.1,4.1,4,2.1))

> hist(d$asc0, main="Approval (Pre) in History", breaks=seq(0.5, 5.5, 1), freq=FALSE, labels=TRUE, col="blue", xaxt='n', xlab="", cex=0.8)

> axis(side=1, at=c(1, 3, 5), labels=c("Strongly \n disapprove", "Neither approve \n nor disapprove", "Strongly \n approve"))

> dev.off()
RStudioGD 
        2 

> hist(d$asc, main="Historical Approval (Post)", breaks=seq(0.5, 5.5, 1), freq=FALSE, labels=TRUE,  col="blue", xaxt='n', xlab="")

> axis(side=1, at=c(1, 3, 5), labels=c("Strongly disapprove", "Neither approve nor disapprove", "Strongly approve"))

> pdf(paste(fig.loc,"Hist-Approval(Post)-Histogram",".pdf", sep=""), width=w, height=h*1.6)

> par(mar=c(3.1,4.1,3,2.1))

> hist(d$asc, main="Approval (Post) in History", breaks=seq(0.5, 5.5, 1), freq=FALSE, labels=TRUE, col="blue", xaxt='n', xlab="", cex=0.6)

> axis(side=1, at=c(1, 3, 5), labels=c("Strongly \n disapprove", "Neither approve \n nor disapprove", "Strongly \n approve"))

> dev.off()
RStudioGD 
        2 

> # #par(oma=c(0,0,-2,0))
> # title("Historical Approval (Post)", line = 2)
> 
> #How many respondents got ADIZp
> sum(d$ADIZp==1 & d$att==1)
[1] 881

> sum(d$ADIZ==1 & d$att==1)
[1] 1758

> table(d$ADIZ, d$ADIZp)
   
       0    1
  0 3410    0
  1 1262 1278

> sum(d$pro==1 & d$att==1)
[1] 858

> pro.results$var <- as.factor(pro.results$var)

> pro.results$dv <- as.factor(pro.results$dv)

> pro.results$cov <- as.factor(pro.results$cov)

> pro.results$est <- as.numeric(pro.results$est)

> pro.results$se <- as.numeric(pro.results$se)

> pro.results$t <- as.numeric(pro.results$t)

> pro.results$p <- as.numeric(pro.results$p)

> pro.results$names2 <- paste(pro.results$var," \n ", pro.results$cov, sep="")

> pro.results
                    est         se          t          p   var       dv           cov                 names2
pro.results -0.05649917 0.03467307 -1.6294828 0.10326880  EP-3 approval no covariates  EP-3 \n no covariates
            -0.03517927 0.04231808 -0.8313057 0.40583744 ADIZp approval no covariates ADIZp \n no covariates
3           -0.05417796 0.03123547 -1.7345009 0.08288609  EP-3 approval    covariates     EP-3 \n covariates
4           -0.04936984 0.03812182 -1.2950547 0.19535653 ADIZp approval    covariates    ADIZp \n covariates
5            1.49055481 0.80948292  1.8413666 0.06562394  EP-3  resolve no covariates  EP-3 \n no covariates
6            1.50050092 0.98718069  1.5199861 0.12857437 ADIZp  resolve no covariates ADIZp \n no covariates
7            0.68126011 0.53321579  1.2776443 0.20143124  EP-3  resolve    covariates     EP-3 \n covariates
8            0.68363434 0.65023350  1.0513675 0.29313813 ADIZp  resolve    covariates    ADIZp \n covariates

> statement.results <- as.data.frame(statement.results)

> statement.results$names <- as.factor(row.names(statement.results))

> statement.results$names2 <- as.factor(rep(c(rep("No Covariates",4),rep("Covariates",4)),2))

> statement.results$names3 <- as.factor(rep(c("Vague Threat", "Nationalist History", "Biding Time", "Cost of War"),2))

> statement.results$names3 <- factor(statement.results$names3, 
+                                    levels=rev(c("Vague Threat", "Nationalist History", "Cost of War", "Biding Time")))

> head(statement.results)
                      Estimate Std. Error   t value   Pr(>|t|)               names        names2              names3
ADIZ, approval      0.05135942 0.03620402 1.4186109 0.15606980      ADIZ, approval No Covariates        Vague Threat
his, approval       0.06902725 0.03425302 2.0152168 0.04393107       his, approval No Covariates Nationalist History
eli.f, approval     0.05999469 0.05418289 1.1072628 0.26822934     eli.f, approval No Covariates         Biding Time
eli.c, approval     0.01414785 0.05547580 0.2550274 0.79871156     eli.c, approval No Covariates         Cost of War
ADIZ, approval, cov 0.06646860 0.03262682 2.0372379 0.04167469 ADIZ, approval, cov    Covariates        Vague Threat
his, approval, cov  0.07167621 0.03088243 2.3209381 0.02032710  his, approval, cov    Covariates Nationalist History

> # factor(x$name, levels = x$name[order(x$val)])
> 
> statement.results$dv <- c(rep(0,8), rep(1,8))

> names(statement.results) <- c("est", "se", "t", "p", "names", "names2","names3", "dv")

> statement.results$names4 <- paste(statement.results$names3, "\n ", statement.results$names2, sep="")

> statement.results
                             est         se           t          p                names        names2              names3 dv
ADIZ, approval        0.05135942 0.03620402  1.41861095 0.15606980       ADIZ, approval No Covariates        Vague Threat  0
his, approval         0.06902725 0.03425302  2.01521676 0.04393107        his, approval No Covariates Nationalist History  0
eli.f, approval       0.05999469 0.05418289  1.10726275 0.26822934      eli.f, approval No Covariates         Biding Time  0
eli.c, approval       0.01414785 0.05547580  0.25502744 0.79871156      eli.c, approval No Covariates         Cost of War  0
ADIZ, approval, cov   0.06646860 0.03262682  2.03723789 0.04167469  ADIZ, approval, cov    Covariates        Vague Threat  0
his, approval, cov    0.07167621 0.03088243  2.32093806 0.02032710   his, approval, cov    Covariates Nationalist History  0
eli.f, approval, cov  0.08783447 0.04881616  1.79929072 0.07202831 eli.f, approval, cov    Covariates         Biding Time  0
eli.c, approval, cov  0.01976154 0.04998947  0.39531402 0.69262671 eli.c, approval, cov    Covariates         Cost of War  0
ADIZ, resolve        -0.50295740 0.84703636 -0.59378490 0.55268144        ADIZ, resolve No Covariates        Vague Threat  1
his, resolve         -0.94916878 0.80005659 -1.18637706 0.23552683         his, resolve No Covariates Nationalist History  1
eli.f, resolve       -0.76095065 1.26122612 -0.60334197 0.54630714       eli.f, resolve No Covariates         Biding Time  1
eli.c, resolve       -0.86972955 1.29806151 -0.67002183 0.50287319       eli.c, resolve No Covariates         Cost of War  1
ADIZ, resolve, cov   -1.00970639 0.55803145 -1.80940767 0.07044473   ADIZ, resolve, cov    Covariates        Vague Threat  1
his, resolve, cov    -0.23212353 0.52756088 -0.43999383 0.65995965    his, resolve, cov    Covariates Nationalist History  1
eli.f, resolve, cov   0.04114682 0.83063093  0.04953682 0.96049337  eli.f, resolve, cov    Covariates         Biding Time  1
eli.c, resolve, cov   0.66982883 0.85542510  0.78303621 0.43364109  eli.c, resolve, cov    Covariates         Cost of War  1
                                                  names4
ADIZ, approval              Vague Threat\n No Covariates
his, approval        Nationalist History\n No Covariates
eli.f, approval              Biding Time\n No Covariates
eli.c, approval              Cost of War\n No Covariates
ADIZ, approval, cov            Vague Threat\n Covariates
his, approval, cov      Nationalist History\n Covariates
eli.f, approval, cov            Biding Time\n Covariates
eli.c, approval, cov            Cost of War\n Covariates
ADIZ, resolve               Vague Threat\n No Covariates
his, resolve         Nationalist History\n No Covariates
eli.f, resolve               Biding Time\n No Covariates
eli.c, resolve               Cost of War\n No Covariates
ADIZ, resolve, cov             Vague Threat\n Covariates
his, resolve, cov       Nationalist History\n Covariates
eli.f, resolve, cov             Biding Time\n Covariates
eli.c, resolve, cov             Cost of War\n Covariates

> pvalues <- -1

> r <- !is.na(d$asc) & rb 

> n <- length(d[r,1])

> #statements, statements cov, mob, mob covs, adiz, adiz covs, 
> 
> type <- "Hist"

> r2 <- statement.results$dv==0 & statement.results$names3=="Vague Threat"   #& statement.results$names2=="No Covariates" 

> dv.n <- "Approval"

> title.w <- "Effect on"

> covs <- ""

> main.t <- "1StateA-"

> ymin.n <- -1; ymax.n <- .4

> source("./MASTER CODE/6-coef-figuresb.R", echo=TRUE)

> #6-coef-figuresb.r
> 
> 
> sign.n <- 1

> hjust.n <- -0.2

> pvalues <- round((statement.results[r2,]$p),3)

> pvalues.labels <- paste(rep("p=",length(pvalues)), pvalues, sep="")

> title.f <- paste(main.t, type, "-", dv.n, sep="")

> title.fig <- paste(title.w, " ", dv.n, ", ", type, sep="")

> if (covs!=""){
+   title.f <- paste(main.t, type, "-", dv.n, "-", covs, sep="")
+   title.fig <- paste(title.w, " ", dv.n, ", ", type, ", ", covs, s .... [TRUNCATED] 

> ggplot(statement.results[r2,]) +    
+   geom_linerange(aes(x = names4, ymin = sign.n * est - outerCI.n*se, ymax = sign.n * est + outerCI.n*se), lwd .... [TRUNCATED] 

> p <- recordPlot()

> pdf(paste(fig.loc,title.f,".pdf", sep=""), width=9, height=3)

> p

> dev.off()
RStudioGD 
        2 

> pvalues <- -1

> type <- "Hist"

> r2 <- statement.results$dv==0 & statement.results$names3!="Vague Threat"   #& statement.results$names2=="No Covariates" 

> dv.n <- "Approval"

> title.w <- "Effect on"

> covs <- ""

> main.t <- "1StateB-"

> ymin.n <- -1; ymax.n <- .4

> source("./MASTER CODE/6-coef-figuresb.R", echo=TRUE)

> #6-coef-figuresb.r
> 
> 
> sign.n <- 1

> hjust.n <- -0.2

> pvalues <- round((statement.results[r2,]$p),3)

> pvalues.labels <- paste(rep("p=",length(pvalues)), pvalues, sep="")

> title.f <- paste(main.t, type, "-", dv.n, sep="")

> title.fig <- paste(title.w, " ", dv.n, ", ", type, sep="")

> if (covs!=""){
+   title.f <- paste(main.t, type, "-", dv.n, "-", covs, sep="")
+   title.fig <- paste(title.w, " ", dv.n, ", ", type, ", ", covs, s .... [TRUNCATED] 

> ggplot(statement.results[r2,]) +    
+   geom_linerange(aes(x = names4, ymin = sign.n * est - outerCI.n*se, ymax = sign.n * est + outerCI.n*se), lwd .... [TRUNCATED] 

> p <- recordPlot()

> pdf(paste(fig.loc,title.f,".pdf", sep=""), width=9, height=3)

> p

> dev.off()
RStudioGD 
        2 

> pvalues <- -1

> type <- "Hist"

> r2 <- statement.results$dv==0 & statement.results$names2=="Covariates"

> dv.n <- "Approval"

> title.w <- "Effect on"

> covs <- "Covs"

> main.t <- "1State-"

> ymin.n <- -1; ymax.n <- .4

> source("./MASTER CODE/6-coef-figuresb.R", echo=TRUE)

> #6-coef-figuresb.r
> 
> 
> sign.n <- 1

> hjust.n <- -0.2

> pvalues <- round((statement.results[r2,]$p),3)

> pvalues.labels <- paste(rep("p=",length(pvalues)), pvalues, sep="")

> title.f <- paste(main.t, type, "-", dv.n, sep="")

> title.fig <- paste(title.w, " ", dv.n, ", ", type, sep="")

> if (covs!=""){
+   title.f <- paste(main.t, type, "-", dv.n, "-", covs, sep="")
+   title.fig <- paste(title.w, " ", dv.n, ", ", type, ", ", covs, s .... [TRUNCATED] 

> ggplot(statement.results[r2,]) +    
+   geom_linerange(aes(x = names4, ymin = sign.n * est - outerCI.n*se, ymax = sign.n * est + outerCI.n*se), lwd .... [TRUNCATED] 

> p <- recordPlot()

> pdf(paste(fig.loc,title.f,".pdf", sep=""), width=9, height=3)

> p

> dev.off()
RStudioGD 
        2 

> pvalues <- -1

> #Demographics 
> 
> 
> rm(list = ls())

> setwd("~/Dropbox/Dafoe-Weiss-Empirics/18-04-12-ISQ-Raluca/Code/Replication File ISQ/")

> load("./MASTER CODE/real.Rdata")

> #Figures location
> fig.loc <- "../../Figures/"

> w <- 5

> h <- 4

> innerCI.n <- 1.64 #set so exclusion of CI implies two sided significance of 0.1

> outerCI.n <- 1.96 #set so exclusion of CI implies two sided significance of 0.05

> shadowtext <- function(x, y=NULL, labels, col='white', bg='black',
+                        theta= seq(pi/4, 2*pi, length.out=8), r=0.1, ... ) {
+   
+   xy <- xy.coords(x,y)
+   xo <- r*strwidth('A')
+   yo <- r*strheight('A')
+   for (i in theta) {
+     text( xy$x + cos(i)*xo, xy$y + sin(i)*yo, 
+           labels, col=bg, ... )
+   }
+   text(xy$x, xy$y, labels, col=col, ... )
+ }

> stat_sum_df <- function(fun, geom="crossbar", ...) {
+   stat_summary(fun.data=fun, colour="red", geom=geom, width=0.2, ...)
+ }

> pdf(paste(fig.loc,"Hist-Approval(Pre)-Hist",".pdf", sep=""), width=w, height=h)

> par(mar=c(3.1,4.1,4,2.1))

> hist(d$asc0, main="Approval (Pre) in Hist", breaks=seq(0.5, 5.5, 1), freq=FALSE, labels=TRUE, col="blue", xaxt='n', xlab="", ylim=c(0,0.5))

> axis(side=1, at=c(1, 3, 5), labels=c("Strongly \n disapprove", "Neither approve \n nor disapprove", "Strongly \n approve"))

> dev.off()
RStudioGD 
        2 

> plot(density(d$rac0[d$pre.questions==1]), main="Resolve")

> #Second, we would expect to see high levels of agreement with the question $cc$: To what extent do you agree or disagree with the following statement: ``Criticism of government policy is unhelpful.'' To the extent that we do not see see overwhelming high levels of $cc$, we gain confidence that respondents are not self-censoring. 
> hist(as.numeric(as.character(d$cc)))

> hist(as.numeric(as.character(d$ccr)))

> summary(lm(d$ccc ~ d$cc.or)) #A small proportion of the answer can be attributed to the order of the answer options. 

Call:
lm(formula = d$ccc ~ d$cc.or)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.1460 -0.9008  0.0992  0.8540  2.0992 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  3.14600    0.02194 143.386  < 2e-16 ***
d$cc.or     -0.24520    0.03107  -7.893 3.72e-15 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.014 on 4262 degrees of freedom
  (1686 observations deleted due to missingness)
Multiple R-squared:  0.01441,	Adjusted R-squared:  0.01418 
F-statistic:  62.3 on 1 and 4262 DF,  p-value: 3.722e-15


> pdf(paste(fig.loc,"Hist-Criticism",".pdf", sep=""), width=w, height=h)

> par(mar=c(3.1,4.1,4,2.1))

> hist(d$ccc, main="Criticism of the Government is Unhelpful", breaks=seq(0.5, 5.5, 1), freq=FALSE, labels=TRUE, col="blue", xaxt='n', xlab="", ylim=c(0,0.5))

> axis(side=1, at=c(1, 3, 5), labels=c("Strongly \n disagree", "Neither disagree \n nor agree", "Strongly \n agree"))

> dev.off()
RStudioGD 
        2 

> #Pre Approval doesn't correlate with support of government criticism
> summary(lm(d$ccc ~ d$asc0)) 

Call:
lm(formula = d$ccc ~ d$asc0)

Residuals:
     Min       1Q   Median       3Q      Max 
-2.05900 -0.98928 -0.00671  0.99329  2.01072 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  3.07644    0.07255  42.403   <2e-16 ***
d$asc0      -0.01743    0.01812  -0.962    0.336    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.015 on 3619 degrees of freedom
  (2329 observations deleted due to missingness)
Multiple R-squared:  0.0002556,	Adjusted R-squared:  -2.064e-05 
F-statistic: 0.9253 on 1 and 3619 DF,  p-value: 0.3362


> #Concern for national honor slightly negatively associated with support of government criticism
> summary(lm(d$ccc ~ d$na1.v)) 

Call:
lm(formula = d$ccc ~ d$na1.v)

Residuals:
     Min       1Q   Median       3Q      Max 
-2.13635 -0.97397 -0.00645  0.97731  2.02603 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  3.136354   0.062696  50.024   <2e-16 ***
d$na1.v     -0.016238   0.007673  -2.116   0.0344 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.015 on 3619 degrees of freedom
  (2329 observations deleted due to missingness)
Multiple R-squared:  0.001236,	Adjusted R-squared:  0.00096 
F-statistic: 4.479 on 1 and 3619 DF,  p-value: 0.03439


> ## Number of observations who passed attention checks
> sum(d$att)
[1] 3717

> sum((d$ctri <= 2) %in% TRUE)
[1] 0

> hist(d$age)

> age <- d$age

> #65% male of attentives; 75% male of initial
> prop.table(table(d$gender))

        0         1 
0.3721008 0.6278992 

> d$educ.f <- as.factor(d$educ)

> d$educ.f[d$educ.m==1] <- NA

> library(plyr)

> # 01 No formal education 02 Elementary school 03 Middle school
> # 04 High school
> # 05 College 06 Masters 07 Doctoral
> d$educ.f <- revalue(d$educ.f, c("1"="No Educ", "2"="Elementary", "3"="Middle", "4"="High School", "5"="College", "6"="MA", "7"="Doctoral"))

> prop.q.data <- table(d$educ.f)/sum(!is.na(d$educ.f))

> prop.ci.data <- c(0,.119,.36,.312,.209, 0, 0)

> prop.data <- rbind(prop.q.data, prop.ci.data)

> row.names(prop.data) <- c("Our Sample", "Chinese \nInternet Users")

> pdf(paste(fig.loc,"hist_education.pdf", sep=""), width=8, height=6)

> barplot(prop.q.data, main="Education", cex.names=0.7)

> dev.off()
RStudioGD 
        2 

> pdf(paste(fig.loc,"hist_education_comp.pdf", sep=""), width=8, height=6)

> barplot(prop.data, main="Education", cex.names=0.7, beside=TRUE, legend=rownames(prop.data))

> dev.off()
RStudioGD 
        2 

> p <- recordPlot()

> p

> d$educ.f2 <- d$educ.f

> d$educ.f2[d$educ.f=="MA" | d$educ.f=="Doctoral"] <- "MA"

> d$educ.f2[d$educ.f=="No Educ"] <- "Elementary"

> d$educ.f2 <- revalue(d$educ.f2, c("MA"="Graduate"))

> d$educ.f2 <- droplevels(d$educ.f2, c("No Educ","Doctoral"))

> prop.q.data <- table(d$educ.f2)/sum(!is.na(d$educ.f2))

> #Primary, Junior, Senior, 2 year college + 4 year college, Graduate
> prop.huang.data <- c(0.006, 0.038, 0.103, 0.307 + 0.507, 0.04)

> prop.data <- rbind(prop.q.data, prop.huang.data)

> row.names(prop.data) <- c("Our Sample", "Huang Data")

> barplot(prop.data, main="Education", cex.names=0.7, beside=TRUE, legend.text=rownames(prop.data), args.legend=list(x=6, y=0.8))

> pdf(paste(fig.loc,"hist_education_comp_Huang.pdf", sep=""), width=8, height=6)

> barplot(prop.data, main="Education", cex.names=0.7, beside=TRUE, legend.text=rownames(prop.data), args.legend=list(x=6, y=0.8))

> dev.off()
RStudioGD 
        2 

> party <- as.factor(d$d13)

> party <- revalue(party, c("1"="Belong now", "2"="Previously", "3"="Have Never", "99"="No Answer"))

> party <- factor(party[party!=""])

> pdf(paste(fig.loc,"hist_CP.pdf", sep=""), width=8, height=6)

> barplot(table(party)/sum(!is.na(party)), main="CP Membership")

> dev.off()
RStudioGD 
        2 

> #### Age analysis ####
> r <- d$age.m==0

> pdf(paste(fig.loc,"hist_age_den.pdf", sep=""), width=8, height=6)

> den <- density(age[r]) # returns the density data 

> plot(den, col="cyan3", main="Age of Respondents") 

> polygon(den, col="cyan3", border="black")

> rug(jitter(age, amount=0.5))

> dev.off()
RStudioGD 
        2 

> ##Actual Chinese age distribution based on 2010 census, eye-ball estimated from this figure:
> #http://en.wikipedia.org/wiki/File:China_Sex_By_Age_2010_census.png
> 
> #From: https://www.cia.gov/library/publications/the-world-factbook/fields/2010.html
> #0-14 years: 17.1% (male 124,340,516/female 107,287,324) 
> #15-24 years: 14.7% (male 105,763,058/female 93,903,845) 
> #25-54 years: 47.2% (male 327,130,324/female 313,029,536) 
> #55-64 years: 11.3% (male 77,751,100/female 75,737,968) 
> #65 years and over: 9.6% (male 62,646,075/female 68,102,830) (2014 est.)
> 
> 
> age.b <- c(0,15,25,55,65,80)

> #Data for China
> age.c <- rep(NA, 1000)

> i <- 1

> j <- i+170

> age.c[i:j] <- runif(171,0,15)  

> i <- j+1

> j <- i+146

> age.c[i:j] <- runif(147,15,25)

> i <- j+1

> j <- i+471

> age.c[i:j] <- runif(472,25,55)

> i <- j+1

> j <- i+113

> age.c[i:j] <- runif(114,55,65)

> i <- j+1

> j <- i+095

> age.c[i:j] <- runif(96,65,80)

> den.c <-density(age.c) # returns the density data 

> #Data for Chinese internet users, from 33rd Statistical Report of Internet Development in China (CNNIC, January 2014).
> #from https://www.google.com/url?sa=t&rct=j&q=&esrc=s&source=web&cd=1&cad=rja&uact=8&ved=0ahUKEwjjrvvUvInNAhVLdz4KHa6VAJsQFggdMAA&url=http%3A%2F%2Fwww1.cnnic.cn%2FIDR%2FReportDownloads%2F201404%2FU020140417607531610855.pdf&usg=AFQjCNEuKJbXON9y_nkjuvfjcjNS445cfg&sig2=gEIJ1w00TJ9tCJCFWmfVmQ
> age.categories.c.i <- c(0, 10,19,29,39,49,59,80)

> age.prop.c.i <- c(2.1,24.5, 30.7, 23.4, 12, 5.2, 2.1)

> ld <- 1000

> age.ci <- NA

> for (i in 1:length(age.prop.c.i)){
+   age.ci <- c(age.ci, runif(round(ld*(age.prop.c.i[i]/100)), age.categories.c.i[i], age.categories.c.i[i+1]))
+ }

> hist(age.ci)

> pdf(paste(fig.loc,"hist_hist_age.pdf", sep=""), width=13, height=8)

> hist(age.c, breaks=c(0,15,25,55,65,80), col=rgb(1,0,0,.5), ylim=c(0,0.03), 
+      xlab="Age", main="Histogram of Age in China (light red) and Qualtrics (green)") 

> hist(age[r], breaks=c(min(age[r]),15,25,55,65,max(age[r])), col=rgb(0,1,0,.5), add=T) 

> dev.off()
RStudioGD 
        2 

> pdf(paste(fig.loc,"hist_den_age_comparison.pdf", sep=""), width=8, height=6)

> hist(age.c, breaks=c(0,15,25,55,65,80), col=rgb(1,0,0,.5), ylim=c(0,0.05), 
+      xlab="Age", main="Histogram/Density of Age in \n China, Chinese Internet, and Qualtrics") 

> text(60,0.015,"Chinese", col=rgb(1,0,0,1))

> #plot(den, col="cyan3", main="Age of Respondents") 
> hist(age.ci, breaks=age.categories.c.i, col=rgb(0,0,1,0.5), add=TRUE)

> text(10,.032, "Chinese Internet Users", col=rgb(0,0,1,1))

> polygon(den, col=rgb(0,1,0,.5), border="black")

> text(45,0.03,"Qualtrics", col=rgb(0,1,0,1))

> rug(jitter(age, amount=0.5), col=rgb(0,1,0,0.5))

> dev.off()
RStudioGD 
        2 

> ###Age data from the UN
> un.age <- read.csv("./MASTER CODE/UNdata_age.csv")

> un.age2 <- subset(un.age, Area=="Total" & Sex=="Both Sexes" & (Age %in% seq(0,99,1)))

> un.age3 <- un.age2[, c("Age", "Value")]

> un.age3$Age <- as.numeric(as.character(un.age3$Age))

> dat <- transform(un.age3, density = Value / sum(Value))

> plot(density ~ Age, data = dat, type = "n")

> lines(density ~ Age, data = dat, col = "red")

> pdf(paste(fig.loc,"hist_den_age_comparison2.pdf", sep=""), width=8, height=6)

> plot(c(0,100),c(0,0.05), ylab="Density", xlab="Age", type="n", main="Density of Age in \n China, Chinese Internet, and Qualtrics")

> text(65,0.015, "China's Age Density", col="red")

> hist(age.ci, breaks=age.categories.c.i, col=rgb(0.2,1,0.2,.5), add=TRUE)

> text(10,.032, "Chinese \n Internet Users", col=rgb(0.2,1,0.2))

> polygon(den, col=rgb(0,0,1,0.5), border="black")

> text(55,0.03, "Qualtrics' Age Density", col=rgb(0,0,1))

> lines(density ~ Age, data = dat, col = "red", lwd=4)

> rug(jitter(age, amount=0.5), col=rgb(0,1,0,0.5))

> dev.off()
RStudioGD 
        2 

> #Ethnicity
> var <- d$d4

> var.r <- d$d4r

> out <- rep(NA, length(var))

> order.r <- rep(NA, length(var))

> r <- !is.na(as.numeric(as.character(var)))

> out[r] <- as.numeric(as.character(var))[r]

> order.r[r] <- 0

> r <- !is.na(as.numeric(as.character(var.r)))

> out[r] <- as.numeric(as.character(var.r))[r]

> order.r[r] <- 1

> summary(lm(out ~ order.r))

Call:
lm(formula = out ~ order.r)

Residuals:
   Min     1Q Median     3Q    Max 
-1.098 -1.098 -0.744 -0.744 97.256 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   1.7435     0.2207   7.901 3.67e-15 ***
order.r       0.3541     0.3124   1.133    0.257    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 9.284 on 3531 degrees of freedom
  (2417 observations deleted due to missingness)
Multiple R-squared:  0.0003637,	Adjusted R-squared:  8.057e-05 
F-statistic: 1.285 on 1 and 3531 DF,  p-value: 0.2571


> d$ethnicity <- out

> d$ethnicity.or <- order.r

> d$ethnicity2 <- as.factor(d$ethnicity)

> d$ethnicity2 <- revalue(d$ethnicity2, c("1"="Han", "2"="Minority", "99"="No Answer"))

> num.han <- round(sum(d$ethnicity2=="Han", na.rm=T)/sum(!is.na(d$ethnicity2)),2)

> pdf(paste(fig.loc,"hist_ethnicity.pdf", sep=""), width=8, height=6)

> barplot(table(d$ethnicity2)/sum(!is.na(d$ethnicity2)), main="Ethnicity")

> text(0.7, 0.8, paste(num.han, "%"))

> dev.off()
RStudioGD 
        2 

> # D2. In what province do you live? 
> 
> region.mapping <- read.csv(file="../../../Original_data/regions.csv")

> d$region <- factor(d$d2, levels=region.mapping$Code, labels=region.mapping$Region)

> # d$region <- factor(d$d2, levels=c(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,99), 
> #                labels=c("Anhui", "Beijing", "Chonqing", "Fujian", "Gansu", "Guangdong", "Guangxi", "Hainan", "Hebei", "Heilongjian",
> #                         "Henan", "Hong Kong", "Hubei", "Hunan", "Inner Mongolia", "Jiangus", "Jiangxi", "Jilin", "Liaoning", "Ningxia", "Qinghai",
> #                         "Shaanxi", "Shandong", "Shanghai", "Shanxi", "Sichuan", "Taiwan", "Tianjin", "Tibet", "Xinjian", "Yunnan", "Zhejian", "No Answer"))
> title <-"In what province do you live?"

> reported.region <- sort(table(d$region))/length(d$region)

> write.table(reported.region, file=paste(fig.loc, "reported_region.csv", sep=""))

> pdf(paste(fig.loc,"region_demographics.pdf", sep=""), width=8*1.3, height=6*1.3)

> par(mar=c(5,8,4,2))

> barplot(reported.region, horiz=TRUE, las=1, cex.names=1, main="Reported Region")

> dev.off()
RStudioGD 
        2 

> ### Analyzing Distribution of na1.v and na2.v
> 
> library(plot3D)

> na1.v <- d$na1.v

> na2.v <- d$na2.v

> d$na2.v.o <-  as.factor(d$na2.v.o)

> # 01 No formal education 02 Elementary school 03 Middle school
> # 04 High school
> # 05 College 06 Masters 07 Doctoral
> d$na2.v.o <- revalue(d$na2.v.o, c("1"="Too Much", "2"="About Right", "3"="Too Little", "8"="Don't Know", "9"="Refuse to Answer"))

> prop <- table(d$na2.v.o)/sum(!is.na(d$na2.v.o))

> prop.other <- c(.169, .273, .36, .186, .013)

> prop.data <- rbind(prop, prop.other)

> row.names(prop.data) <- c("Our Sample", "RCCC 2012")

> prop.data
            Too Much About Right Too Little Don't Know Refuse to Answer
Our Sample 0.1016489   0.3414736  0.4758923 0.06720935       0.01377583
RCCC 2012  0.1690000   0.2730000  0.3600000 0.18600000       0.01300000

> barplot(prop.data, main="Does China rely on military strength ... (na2)", cex.names=0.7, beside=TRUE, legend=rownames(prop.data))

> pdf(paste(fig.loc,"na2comp1.pdf", sep=""), width=8, height=6)

> barplot(prop.data, main="Does China rely on military strength ... (na2)", cex.names=0.7, beside=TRUE, legend=rownames(prop.data))

> dev.off()
RStudioGD 
        2 

> prop2 <- as.data.frame(prop/sum(prop[1:3]))[1:3,]

> prop.other2 <- prop.other[1:3]/sum(prop.other[1:3])

> prop.data2 <- cbind(prop2, prop.other2)

> colnames(prop.data2) <- c("", "Our Sample", "RCCC 2012")

> prop.data2 <- t(prop.data2)

> colnames(prop.data2) <- prop.data2[1,]

> prop.data2 <- prop.data2[2:3,]

> class(prop.data2) <- "numeric"

> barplot(prop.data2, main="Does China rely on military strength ... (na2)", cex.names=0.7, beside=TRUE, legend=rownames(prop.data), args.legend=list(x="topleft"))

> pdf(paste(fig.loc,"na2comp2.pdf", sep=""), width=8, height=6)

> barplot(prop.data2, main="Does China rely on military strength ... (na2)", cex.names=0.7, beside=TRUE, legend=rownames(prop.data), args.legend=list(x="topleft"))

> #x=5, y=.55))
> dev.off()
RStudioGD 
        2 

> ##  Create cuts:
> x_c <- cut(na1.v, 11)

> y_c <- cut(na2.v, 3)

> ##  Calculate joint counts at cut levels:
> z <- table(x_c, y_c)

> ##  Plot as a 3D histogram:
> hist3D(z=z, border="black", xlab="\n Support for \n National Honor (na1)", ylab="\n China Should \n Use More Force (na2)", theta=-40, phi=20, main="Distribution of Pre-Scenario Responses")

> p <- recordPlot()

> pdf("../../Figures/na1na2.pdf", width=1.5*w, height=1.5*h)

> p

> dev.off()
RStudioGD 
        2 

> ## Analyzing aso ###
> 
> rm(list=ls())

> d0 <- read.csv("./MASTER CODE/real.csv")

> d1 <- read.csv("../../../Original_data/aso/150515RealHistorys-aso.csv")

> varstokeep <- c("V1", "reputation",	"strength/action","sovereignty", "honor", "biding.time", "peace", "domestic", "interests", "nationalism", "deference", "complex", "nonsense.missing")

> d2 <- d1[, names(d1) %in% varstokeep]

> d0$V1<- d0$V1

> d <- merge(d0,d2, by="V1")

> remove(d0, d1, d2)

> asovars <- varstokeep[-1]

> asovarsindex <- names(d) %in% asovars

> m <- is.na(d[, asovarsindex])

> d[, asovarsindex][m] <- 0

> aso.means <- colMeans(d[, asovarsindex])

> aso.means <- aso.means[1:10]

> aso.means.n <- as.numeric(aso.means)[1:10]

> w <- 5

> h <- 4

> # Fitting Labels 
> par(las=2) # make label text perpendicular to axis

> par(mar=c(5,8,4,2)) # increase y-axis margin.

> bp <- barplot(aso.means, xlab="Proportion ", horiz=TRUE, width=0.7, main="Reasons Given")

> text(x=aso.means.n-0.006, y=bp, round(aso.means, 2),cex=1) 

> p <- recordPlot()

> pdf("../../Figures/aso1.pdf", width=1.5*w, height=1.5*h)

> p

> dev.off()
RStudioGD 
        2 

> #-------------------------------------------------------------------------
> #National Honour (Prequestions) Replication
> #-------------------------------------------------------------------------
> 
> rm(list = ls())

> # Customized Stargazer Function, without asterisk inflation
> stargazerAD <- function(table, title="Default", filename="Figures/default.tex", dep.var.labels=NULL) {
+   a <- stargazer(table, title=title,  out=filename, 
+                  align=TRUE, star.cutoffs = c(0.1, 0.05, 0.01, 0.001), star.char = c("\\dagger", "*", "**", "***"), 
+                  notes="$^{\\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$", notes.append=FALSE, single.row=TRUE)
+   return(a)
+ }

> setwd("~/Dropbox/Dafoe-Weiss-Empirics/18-04-12-ISQ-Raluca/Code/Replication File ISQ/")

> d <- read.csv("./MASTER CODE/realanalysis.csv") 

> #when prequestions asked?
> 
> ggplot(d,aes(x=start.time.n,y=pre.questions))+stat_smooth(se=F) + 
+   geom_vline(xintercept = d$start.time.n[1850])

> df <- d[d$start.time.n<d$start.time.n[1850],]

> m <- lm(asc~pre.questions, data=df)

> stargazerAD(m, title="Prequestions Effect on Approval (Real History)", file="../../Figures/real_preq.tex")

% Table created by stargazer v.5.2.2 by Marek Hlavac, Harvard University. E-mail: hlavac at fas.harvard.edu
% Date and time: Fri, Feb 22, 2019 - 09:46:46
% Requires LaTeX packages: dcolumn 
\begin{table}[!htbp] \centering 
  \caption{Prequestions Effect on Approval (Real History)} 
  \label{} 
\begin{tabular}{@{\extracolsep{5pt}}lD{.}{.}{-3} } 
\\[-1.8ex]\hline 
\hline \\[-1.8ex] 
 & \multicolumn{1}{c}{\textit{Dependent variable:}} \\ 
\cline{2-2} 
\\[-1.8ex] & \multicolumn{1}{c}{asc} \\ 
\hline \\[-1.8ex] 
 pre.questions & -0.193^{***}$ $(0.052) \\ 
  Constant & 3.576^{***}$ $(0.041) \\ 
 \hline \\[-1.8ex] 
Observations & \multicolumn{1}{c}{1,843} \\ 
R$^{2}$ & \multicolumn{1}{c}{0.007} \\ 
Adjusted R$^{2}$ & \multicolumn{1}{c}{0.007} \\ 
Residual Std. Error & \multicolumn{1}{c}{1.076 (df = 1841)} \\ 
F Statistic & \multicolumn{1}{c}{13.910$^{***}$ (df = 1; 1841)} \\ 
\hline 
\hline \\[-1.8ex] 
\textit{Note:}  & \multicolumn{1}{r}{$^{\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$} \\ 
\end{tabular} 
\end{table} 
 [1] ""                                                                                                                  
 [2] "% Table created by stargazer v.5.2.2 by Marek Hlavac, Harvard University. E-mail: hlavac at fas.harvard.edu"       
 [3] "% Date and time: Fri, Feb 22, 2019 - 09:46:46"                                                                     
 [4] "% Requires LaTeX packages: dcolumn "                                                                               
 [5] "\\begin{table}[!htbp] \\centering "                                                                                
 [6] "  \\caption{Prequestions Effect on Approval (Real History)} "                                                      
 [7] "  \\label{} "                                                                                                      
 [8] "\\begin{tabular}{@{\\extracolsep{5pt}}lD{.}{.}{-3} } "                                                             
 [9] "\\\\[-1.8ex]\\hline "                                                                                              
[10] "\\hline \\\\[-1.8ex] "                                                                                             
[11] " & \\multicolumn{1}{c}{\\textit{Dependent variable:}} \\\\ "                                                       
[12] "\\cline{2-2} "                                                                                                     
[13] "\\\\[-1.8ex] & \\multicolumn{1}{c}{asc} \\\\ "                                                                     
[14] "\\hline \\\\[-1.8ex] "                                                                                             
[15] " pre.questions & -0.193^{***}$ $(0.052) \\\\ "                                                                     
[16] "  Constant & 3.576^{***}$ $(0.041) \\\\ "                                                                          
[17] " \\hline \\\\[-1.8ex] "                                                                                            
[18] "Observations & \\multicolumn{1}{c}{1,843} \\\\ "                                                                   
[19] "R$^{2}$ & \\multicolumn{1}{c}{0.007} \\\\ "                                                                        
[20] "Adjusted R$^{2}$ & \\multicolumn{1}{c}{0.007} \\\\ "                                                               
[21] "Residual Std. Error & \\multicolumn{1}{c}{1.076 (df = 1841)} \\\\ "                                                
[22] "F Statistic & \\multicolumn{1}{c}{13.910$^{***}$ (df = 1; 1841)} \\\\ "                                            
[23] "\\hline "                                                                                                          
[24] "\\hline \\\\[-1.8ex] "                                                                                             
[25] "\\textit{Note:}  & \\multicolumn{1}{r}{$^{\\dagger} p<0.1$; $^{*} p<0.05$; $^{**} p<0.01$; $^{***} p<0.001$} \\\\ "
[26] "\\end{tabular} "                                                                                                   
[27] "\\end{table} "                                                                                                     

> #-------------------------------------------------------------------------
> #Raw Results
> #-------------------------------------------------------------------------
> rm(list = ls())

> library(ggplot2)

> library(coefplot)

> library(stargazer)

> library(scales)

> library(dplyr)

> detach("package:plyr", unload=TRUE)

> setwd("~/Dropbox/Dafoe-Weiss-Empirics/18-04-12-ISQ-Raluca/Code/Replication File ISQ/")

> d <- read.csv("./MASTER CODE/realanalysis.csv")

> full_data <- d

> myvars <- names(full_data) %in% c("his", "pro", "ADIZ", "ADIZp", "eli.f", "eli.c", "asc") 

> small_data <- full_data[myvars]

> small_data$ID <- seq.int(nrow(small_data))

> #Constructing Appropriate Comparisons 
> #Treatments are assigned independently, except vague threat and defiance occurs only if vague threat
> # only one of biding time or cost of war can occur 
> 
> # the control is those without biding time and those without cost of war (only one can happen)
> small_data$new_eli.f[small_data$eli.f==1 ] <- 1

> small_data$new_eli.f[small_data$eli.f==0 & small_data$eli.c==0] <- 0

> #similarly for cost of war
> 
> small_data$new_eli.c[small_data$eli.c==1] <- 1

> small_data$new_eli.c[small_data$eli.c==0 & small_data$eli.f==0] <- 0

> #making a dataset for each treatment 
> 
> #history
> histvars <- names(small_data) %in% c("his", "asc")

> hist_data <- small_data[histvars]

> hist_data$his_type[small_data$his==1] <- "Treatment"

> hist_data$his_type[small_data$his==0] <- "Control"

> sum(is.na(hist_data$his)) ## no missing for history treatment
[1] 0

> sum(is.na(hist_data$asc)) ## some missing for asc. we'll take them out for the ggplot
[1] 505

> hist_data <- na.omit(hist_data) ##taking out missing values 

> hist_data2 <- hist_data %>% 
+   group_by(his,asc) %>% 
+   summarise(count=n()) %>% 
+   mutate(perc=count/sum(count))

> hist_data2$overall = "History"

> hist_gg <-ggplot(hist_data2, aes(x=as.factor(asc), y=perc, group=as.factor(his), fill=as.factor(his)))

> hist_gg2 <- hist_gg +
+   geom_bar(stat= "identity", position="dodge") + 
+   scale_y_continuous(labels=percent) + 
+   facet_wrap(~overall)

> hist_gg3 <- hist_gg2 +  
+   ylab("Percentage") + 
+   xlab (" ")  + 
+   theme_bw() + 
+   theme(axis.text.x = element_text(angle = 45, hjust = 1))   + 
+   scale_x_discrete(breaks=c(1,2,3,4,5),labels= c("S Disaprove", "Disaprove", "Neither", "Approve", "S Approve")) + 
+   scale_fill_manual(labels=c("Control", "Treatment"), values=c("goldenrod1", "grey40"), name = "Treatment Status") + 
+   coord_cartesian(ylim = c(0,.5))

> #Vague Threat
> 
> ADIZvars <- names(small_data) %in% c("ADIZ", "asc")

> ADIZ_data <- small_data[ADIZvars]

> ADIZ_data$ADIZ_type[small_data$ADIZ==1] <- "Treatment"

> ADIZ_data$ADIZ_type[small_data$ADIZ==0] <- "Control"

> sum(is.na(ADIZ_data$ADIZ_type)) ## no missing 
[1] 0

> sum(is.na(ADIZ_data$asc)) ## some missing for asc. we'll take them out for the ggplot
[1] 505

> ADIZ_data <- na.omit(ADIZ_data) ##taking out missing values 

> ADIZ_data2 <- ADIZ_data %>% 
+   group_by(ADIZ,asc) %>% 
+   summarise(count=n()) %>% 
+   mutate(perc=count/sum(count))

> ADIZ_data2$overall = "Vague Threat"

> ADIZ_gg <-ggplot(ADIZ_data2, aes(x=as.factor(asc), y=perc, group=as.factor(ADIZ), fill=as.factor(ADIZ)))

> ADIZ_gg2 <- ADIZ_gg +
+   geom_bar(stat="identity", position="dodge") + 
+   scale_y_continuous(labels=percent) + 
+   facet_wrap(~overall)

> ADIZ_gg3 <- ADIZ_gg2 +  
+   ylab(" ") + 
+   xlab (" ")  + 
+   theme_bw() + 
+   theme(axis.text.x = element_text(angle = 45, hjust = 1))   + 
+   scale_x_discrete(breaks=c(1,2,3,4,5),labels= c("S Disaprove", "Disaprove", "Neither", "Approve", "S Approve")) + 
+   scale_fill_manual(labels=c("Control", "Treatment"), values=c("goldenrod1", "grey40"), name = "Treatment Status") + 
+   coord_cartesian(ylim = c(0,.5))

> #Biding Time (eli.f)
> 
> elifvars <- names(small_data) %in% c("new_eli.f", "asc")

> elif_data <- small_data[elifvars]

> elif_data$elif_type[small_data$new_eli.f==1] <- "Treatment"

> elif_data$elif_type[small_data$new_eli.f==0] <- "Control"

> sum(is.na(elif_data$elif_type)) ## no missing 
[1] 401

> sum(is.na(elif_data$asc)) ## some missing for asc. we'll take them out for the ggplot
[1] 505

> elif_data <- na.omit(elif_data) ##taking out missing values 

> elif_data2 <- elif_data %>% 
+   group_by(new_eli.f,asc) %>% 
+   summarise(count=n()) %>% 
+   mutate(perc=count/sum(count))

> elif_data2$overall = "Biding Time"

> elif_gg <-ggplot(elif_data2, aes(x=as.factor(asc), y=perc, group=as.factor(new_eli.f), fill=as.factor(new_eli.f)))

> elif_gg2 <- elif_gg +
+   geom_bar(stat="identity", position="dodge") + 
+   scale_y_continuous(labels=percent) + 
+   facet_wrap(~overall)

> elif_gg3 <- elif_gg2 +  
+   ylab("Percentage") + 
+   xlab ("Government Approval")  + 
+   theme_bw() + 
+   theme(axis.text.x = element_text(angle = 45, hjust = 1))   + 
+   scale_x_discrete(breaks=c(1,2,3,4,5),labels= c("S Disaprove", "Disaprove", "Neither", "Approve", "S Approve")) + 
+   scale_fill_manual(labels=c("Control", "Treatment"), values=c("goldenrod1", "grey40"), name = "Treatment Status") + 
+   coord_cartesian(ylim = c(0,.5))

> #Cost of war (eli.c)
> 
> elicvars <- names(small_data) %in% c("new_eli.c", "asc")

> elic_data <- small_data[elicvars]

> elic_data$elic_type[small_data$new_eli.c==1] <- "Treatment"

> elic_data$elic_type[small_data$new_eli.c==0] <- "Control"

> sum(is.na(elic_data$new_eli.c)) ## no missing 
[1] 423

> sum(is.na(elic_data$asc)) ## some missing for asc. we'll take them out for the ggplot
[1] 505

> elic_data <- na.omit(elic_data) ##taking out missing values 

> elic_data2 <- elic_data %>% 
+   group_by(new_eli.c,asc) %>% 
+   summarise(count=n()) %>% 
+   mutate(perc=count/sum(count))

> elic_data2$overall = "Cost of War"

> elic_gg <-ggplot(elic_data2, aes(x=as.factor(asc), y=perc, group=as.factor(new_eli.c), fill=as.factor(new_eli.c)))

> elic_gg2 <- elic_gg +
+   geom_bar(stat="identity", position="dodge") + 
+   scale_y_continuous(labels=percent) + 
+   facet_wrap(~overall)

> elic_gg3 <- elic_gg2 +  
+   ylab(" ") + 
+   xlab ("Government Approval")  + 
+   theme_bw() + 
+   theme(axis.text.x = element_text(angle = 45, hjust = 1))   + 
+   scale_x_discrete(breaks=c(1,2,3,4,5),labels= c("S Disaprove", "Disaprove", "Neither", "Approve", "S Approve")) + 
+   scale_fill_manual(labels=c("Control", "Treatment"), values=c("goldenrod1", "grey40"), name = "Treatment Status") + 
+   coord_cartesian(ylim = c(0,.5))

> #putting them together
> 
> library(ggpubr)

> ggarrange(hist_gg3, ADIZ_gg3, elif_gg3, elic_gg3, 
+           ncol = 2, nrow = 2,common.legend = TRUE, legend = "right")

> ggsave(file="../../Figures/real_vertical_only4.eps")

> #--------------------------------------------------------------------------------------- 
> #Coefficients Plot from Regression
> 
> #Regression Plot
> 
> treat1 <- lm(asc~his, data=hist_data)

> treat3 <- lm(asc~ADIZ, data=ADIZ_data)

> treat5 <- lm(asc~new_eli.f, data=elif_data)

> treat6 <- lm(asc~new_eli.c, data=elic_data)

> mlplot <- multiplot(treat1, treat3, treat5, treat6, 
+                     predictors = c("his", "ADIZ", "new_eli.f", "new_eli.c"), title="Effect on Approval, Hist", xlab="Change in Approval")

> mlplot + theme_bw() + 
+   scale_y_discrete(label=c("History", "Vague Threat", 
+                            "Biding Time", "Cost of War")) + 
+   theme(legend.position="none") + 
+   scale_color_manual(values=c("black", "black", "black", "black", "black", "black"))

> ggsave(file="../../Figures/real_raw_coefplot_only4.eps")

> #Error Bar Plot
> 
> 
> #History
> 
> gg1 <- ggplot(hist_data, aes(as.factor(his_type), asc, fill=as.factor(his_type)))+ 
+   stat_summary(fun.y=mean, geom="bar", width=1, position= "dodge") +
+   stat_summary(fun.data =mean_se, geom="errorbar",  width=.2, size=1, color="black",
+                position=position_dodge(width = .7)) +
+   ylab(" ") +   xlab("History") + 
+   scale_fill_manual(labels=c("Control", "Treatment"), values=c("goldenrod1", "grey40"), name="Treatment Status") + 
+   coord_cartesian(ylim=c(3.2,3.6)) + 
+   theme(panel.background = element_blank(), 
+         axis.ticks.length=unit(0.3,"cm"),
+         axis.ticks.x = element_blank(),
+         axis.text.x = element_blank(), 
+         axis.line.x = element_blank(),
+         axis.line.y = element_line(colour = "Black", linetype = "solid"))+
+   scale_y_continuous(limits = c(3.2,3.6),  expand = c(0,0), oob=rescale_none) 

> #Vague Threat
> 
> gg3 <- ggplot(ADIZ_data, aes(as.factor(ADIZ_type), asc, fill=as.factor(ADIZ_type)))+ 
+   stat_summary(fun.y=mean, geom="bar", width=1, position= "dodge") +
+   stat_summary(fun.data =mean_se, geom="errorbar",  width=.2, size=1, color="black",
+                position=position_dodge(width = .7)) +
+   ylab(" ") +   xlab("VT") + 
+   scale_fill_manual(labels=c("Control", "Treatment"), values=c("goldenrod1", "grey40"), name="Treatment Status") + 
+   scale_y_continuous(limits = c(3.2,3.6),  expand = c(0,0), oob=rescale_none)  + 
+   coord_cartesian(ylim=c(3.2,3.6)) + 
+   theme(panel.background = element_blank(), 
+         axis.line.y = element_blank(),
+         axis.ticks.length=unit(0.3,"cm"),
+         axis.ticks.x = element_blank(),
+         axis.text.x = element_blank(), 
+         axis.ticks.y =element_blank(), 
+         axis.text.y = element_blank(),
+         axis.line.x = element_blank())

> #Biding Time
> 
> gg5 <- ggplot(elif_data, aes(as.factor(elif_type), asc, fill=as.factor(elif_type)))+ 
+   stat_summary(fun.y=mean, geom="bar", width=1, position= "dodge") +
+   stat_summary(fun.data =mean_se, geom="errorbar",  width=.2, size=1, color="black",
+                position=position_dodge(width = .7)) +
+   ylab(" ") +   xlab("Biding") + 
+   scale_fill_manual(labels=c("Control", "Treatment"), values=c("goldenrod1", "grey40"), name="Treatment Status") + 
+   scale_y_continuous(limits = c(3.2,3.6),  expand = c(0,0), oob=rescale_none)  + 
+   coord_cartesian(ylim=c(3.2,3.6)) + 
+   theme(panel.background = element_blank(), 
+         axis.line.y = element_blank(),
+         axis.ticks.length=unit(0.3,"cm"),
+         axis.ticks.x = element_blank(),
+         axis.text.x = element_blank(), 
+         axis.ticks.y =element_blank(), 
+         axis.text.y = element_blank(),
+         axis.line.x = element_blank())

> #Cost of War 
> 
> gg6 <- ggplot(elic_data, aes(as.factor(elic_type), asc, fill=as.factor(elic_type)))+ 
+   stat_summary(fun.y=mean, geom="bar", width=1, position= "dodge") +
+   stat_summary(fun.data =mean_se, geom="errorbar",  width=.2, size=1, color="black",
+                position=position_dodge(width = .7)) +
+   ylab(" ") +   xlab("Cost of War") + 
+   scale_fill_manual(labels=c("Control", "Treatment"), values=c("goldenrod1", "grey40"), name="Treatment Status") + 
+   scale_y_continuous(limits = c(3.2,3.6),  expand = c(0,0), oob=rescale_none)  + 
+   coord_cartesian(ylim=c(3.2,3.6)) + 
+   theme(panel.background = element_blank(), 
+         axis.line.y = element_blank(),
+         axis.ticks.length=unit(0.3,"cm"),
+         axis.ticks.x = element_blank(),
+         axis.text.x = element_blank(), 
+         axis.ticks.y =element_blank(), 
+         axis.text.y = element_blank(),
+         axis.line.x = element_blank())

> ggarrange(gg1, gg3,  gg5, gg6, 
+           ncol = 4, nrow = 1,common.legend = TRUE, legend = "right")

> ggsave(file="../../Figures/real_subst_only4.eps")

> #End. 
